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WALL TEMPERATURE DISTRIBUTION CALCULATION FOR 
A ROCKET NOZZLE CONTOUR 

SUMMARY 


A concept is presented which allows the calculation of the temperatures 
along a thrust chamber nozzle contour on the hot gas and on the coolant flow- 
side. Also considered is a regenerative coolant flowing in the opposite or the 
same direction to the chamber reaction products. Coupling of the boundary 
layer equations for the hot gas-side with the regenerative cooling equations pro- 
vides the results. Since the new analytical model has been integrated into the 
JANNAF Turbulent Boundary Layer computer program, the thrust degradation, 
due to the viscous effects close to the wall, is simultaneously obtained. The 
calculation is started with approximated temperature distributions for the hot 
gas-side wall and the coolant flow. Iterations within the computer program are 
executed until the heat transfer rates from the boundary layer to the wall and 
from the wall to the coolant are equal. Kinetic inviscid flow conditions for 
the boundary layer edge are considered by means of table inputs representing 
the variation of appropriate parameters. Since the chamber wall thickness 
and the coolant flow channel geometry are part of the analysis, optimization 
studies can be performed for these parameters by consecutive computer runs. 

A sample calculation, utilizing the new concept for a small area ratio high 
pressure thrust chamber, is included. 


INTRODUCTION 


The calculation of various turbulent boundary layer thicknesses in the 
thrust chamber and the temperatures of the gas-side wall, the regenerative 
coolant-side wall, and the coolant fluid along the thrust chamber contour is 
simultaneously made, by considering the heat exchange between the combustion 
product flow in the thrust chamber and the coolant flow in the cooling jacket. 

The Turbulent Boundary Layer Computer Program TBL-I [1] has been modi- 
fied to carry out the calculation by using a new concept in which the boundary 
layer equations are coupled with regenerative cooling equations. 

The steady-state conditions that are considered require the temperatures 
of the combustion products, the chamber walls, and also the heat flux through 



the walls to remain constant at any point in time. It is assumed that heat 
transfer occurs only by convection and conduction from the hot combustion 
products to the thrust chamber wall, neglecting the radiation. However, 
inclusion of radiation is not difficult, if the emissivity of the combustion pro- 
ducts and the Stefan- Boltzmann constant can be accurately determined; since 
the total specific heat flux from the hot gas into the chamber wall is composed 

of the convective q and the radiant q heat flux. 

c r 


The coolant fluid in this analysis flows through the tubes or channels in the 
opposite or same direction to the combustion products, receiving the heat by 
convection and conduction. The heat exchange takes place simultaneously in 
many small sections which have an arbitrary length along the contour of the 
thrust chamber, accounting for the gas phase turbulent combustion product 
flow, the temperature of the thrust chamber wall material, and the temperature 
of regenerative coolant flow. The temperature distributions obtained from 
the first iteration are internally used as initial values for the second iteration. 
Iterations are performed until convergence is obtained. Since the influence of 
the coolant transport properties on the resulting temperatures is quite signifi- 
cant, it is important to use pertinent values especially in the supercritical 
region of the coolant fluid. The empirical relationship of the heat transfer co- 
efficient for computing the heat exchange with the coolant flow significantly af- 
fects the results as well as the Stanton number of the combustion products [1], 

The methods to calculate the various turbulent boundary layer thick- 
nesses in the thrust chamber are explained in detail in the documentations of 
TBL [2] and TBL-I [l] , and only the fundamental equations and concepts of 
the calculations improving the latter TBL-I computer program are outlined 
in this report. The concept is demonstrated for a regenerative coolant flowing 
in opposite direction to the combustion products. The alternate equations for 
the coolant flowing in the same direction are explained in the section entitled 
Same Direction Coolant Flow. 
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FUNDAMENTAL EQUATIONS FOR THE BOUNDARY LAYER 


The integral momentum and energy equations in axisymmetric form 
[2,3] for compressible turbulent boundary layer flow are: 



where the displacement thickness 6 , momentum thickness 9 and energy 

thickness <p are identified as follows: 
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The skin friction coefficient is defined as 
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C = 
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f p U 2 


( 6 ) 


and has a form of the Blasius relation [ 1] 
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with the following Reynolds number based upon the momentum thickness 
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The Stanton number 


C_ _ = 




H p U (H - H ) 

co co aw w 


is calculated using the formula [l] 


( 9 ) 


4 




Velocity and enthalpy profiles across the boundary layer are assumed to follow 
the relationships: 


for y — 6 


for y > 6 



( 11 ) 


( 12 ) 
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The definition of enthalpy is 
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h 0 = 
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The adiabatic wall enthalpy H is defined as 

aw 


aw 

hT 


H 

00 



(18) 


The density p within the boundary layer is obtained from the perfect gas equa- 
tion, assuming that the pressure is constant across the boundary layer: 


T 

00 

T 


(19) 


where the temperature T is calculated via the velocity and enthalpy distribu- 
tions, equations (11), (12), (13), (14), and (16). The boundary layer calcu- 
lations use the Runge-Kutta Gill solution method for given parameters at the 
boundary layer edge such as x, r, M , P , T , U , and 3ft. The only 

unknown parameter is the wall temperature T in equation (17) . 


EQUATIONS FOR THE REGENERATIVE COOLING CYCLE 

As shown in Figure 1, the coolant flows in an opposite direction to the 
combustion products of the thrust chamber. The regenerative fluid enters 


6 


downstream with a lower temperature and a higher pressure than at the injector 

head, since heat is continuously transferred from the combustion products to 

the coolant through the chamber walls. T denotes the gas- side wall 

wg 

temperature, the coolant-side wall temperature, and T the coolant 

i 

bulk temperature at an arbitrary station x with x = 0 at the throat (Figs. 

1 and 2) . We consider the case in which the heat is transferred only by con- 
vection from the hot combustion products to the chamber wall and that the 
direction of heat flow is normal to it. Since steady-state conditions are treated, 
the temperatures of the combustion gas and wall and the specific heat flux 
through the walls remain constant with time at any given point. 

The five fundamental equations representing the cooling cycle, includ- 
ing an empirical relation for the heat transfer coefficient of the coolant, are 
as follows: 


1. Specific heat transfer rate on the gas-side, 


q 

wl 




( 20 ) 


where h is the heat transfer coefficient in a gas and T is the adiabatic 
g aw 

wall temperature. 

2. Heat transfer coefficient in a gas related to the Stanton number which 
is calculated in TBL-I, 


h 
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= P„ 


U 

00 



( 21 ) 


where 

is the- free stream density, 
is the free stream velocity, 
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C is the Stanton number, 
H 


H is the adiabatic wall enthalpy, 
aw 

H is the wall enthalpy, 
w 

T is the adiabatic wall temperature, 
aw 

and 

|T wg J. is the input wall temperature or calculated wall temperature. 

3. Specific heat transfer rate through the wall by conduction, 


q 

w2 



T - T 
wg w i 


t 


( 22 ) 


where \ is the thermal conductivity of the wall material and t is the wall 
w 

thickness. 

4. Specific heat transfer rate into the coolant, 



where h^ is the heat transfer coefficient for the coolant. 

5. Empirical relation of the heat transfer coefficient for the hydrogen 
coolant flow [4] is a modified Colburn equation. For any other coolant flow, 
a similar relationship must be utilized including the effects of curvature, 
associated turbulence, and surface roughness of the tubes represented by the 
enhancement factor rj . The accuracy of the enhancement factor significantly 

ht 
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affects the heat transfer calculation and the resulting wall temperatures. 
Since this effect is coupled with the cooling fluid heat transfer coefficient, it 
is evident that the physical property information must be very precise. 


\ / r j l \ 0. 55 

V °' 025 i~ R e * 8 P °r * 4 (t^I 

tube i SL \ w 


(24) 


The above equation is valid for temperature ratios T /T between 1.44 to 

Jt 

9.2, where the Reynolds number and the Prandtl number of the coolant are 
defined as follows: 


Reynolds number, R 


P i U i D tube 


"i M, 


(25) 
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Prandtl number, P = — - — £ — 
r i 


(26) 


Mass flow density, p = P^( x ) u _f( x ) 


(27) 


Equivalent tube diameter, = 2 1 / 2 


(28) 


Coolant bulk viscosity, p^ = p^ , Pressure ^ . (29) 
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Coolant bulk specific heat, 

Coolant bulk thermal conductivity, \ = 


, Pressure j . (30) 

, Pressure^ . (31) 


For steady- state conditions, the heat flux through all three realms must 
be constant. 


V = Vi ' \ 3 * V = constant 


(32) 


Unknowns in equations (20) through (24) are 




T , 
wg 


In equation (21) h is independently calculated when T 
bining equations (20), (22), (23), and (32) results in 


T and T . 
w t 

i 

is given. Com- 



and 


T 

wg 


\ 

, m w _ 

h T + —7 T 
g aw t w. 
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(33) 


(34) 


Derivations of the above equations are shown in Appendix A. Thus, the solu- 
tion can be obtained by considering equations (20) , (21), (24), (33), and (34), 

(Table l). The flow chart to compute T (x) , T (x) , T. (x) , and 

w { wgc v ' ic 
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(x) is shown in Figure 3 where the subscript c denotes the internally 
calculated temperature. 

At the beginning of the calculation, the coolant bulk temperature distri- 
bution is approximated. The coolant-side wall temperature at an axial 

distance x is obtained according to iterations in statement (2) of Figure 3a. 
The gas-side wall temperature is calculated from equation (34); how- 

ever, to differentiate between the input table values of T , the subscript 

wg 

c is added in statement (3) of Figure 3a. The term q ' , which differs from 

w 

the q^ output of TBL-I 

plete. In statement (5) of Figure 3a the coolant temperature is calculated by 
using its previous iteration values. The derivation of the equation in statement 
(5) is shown in the section entitled Internally Calculated Coolant Bulk Tempera- 
ture. 


, should coincide with q^ after the iteration is cora- 


After obtaining T^^x) and T^(x) at each table point of x, the 

values of T„(x) and T (x) to be input for a successive iteration are deter- 
V ’ wg v ’ 

mined as follows: 



T ic<*>+ P 1 /*)], 


(35) 


and 



(36) 


In repeating the preceding calculation, we obtain values of LV x) J 3 and 
^(x)] 3 . This operation is applied until a desired convergence outlined 
in a later section is achieved. 
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INTERNALLY CALCULATED COOLANT BULK TEMPERATURE 


For simplicity, assume that the inner wall of the thrust chamber con- 
sists of a single wall and not of tubes. Let us consider an arbitrary section i 
in Figure 4 and calculate the coolant temperature at x. , which is the distance 

along the nozzle axis. Section i contains the surface area between B and D, 
as shown in Figure 4, 

x = x - A x. , which is x, , (37) 

i - 1 l ll A 


and 


x = x + A x. , which is x^ 
i + 1 i i2 E 


(38) 


where the step sizes A x., and A x. are arbitrary. 

ll i2 


A x. 


The inlet temperature of the coolant at section i is lx. + 




A x 


and the outlet temperature is lx. - 
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. The heat transfer rate through 


the cylindrical surface area of section i between B and D is 


Q (x.) = 

w v 1 


r(x.) qj (x.) Ax. 


COS G! (x.) 
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(39) 


where 


Ax. 
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Ax. 


il 


+ Ax 

i2 
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(40) 
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and a (x.) is the angle between the chamber wall and the nozzle axis at x. . 
The wall radius is r(x^) and q^'( x^ ) is the specific heat transfer rate as 
shown in equation (32). The outlet temperature of the coolant at 

AXii 

x = x in section i is calculated bv 

i 2 J 



T.(x.)| •+ [t,(x. + Ax.j"| . Q (x.) 

£ i J 1 L l i i2 ; J i w v i ; 


+ m C (x.) ’ (41) 
i pi v i' 


where jT^(x.)J ^ and |t^(x. + A x.^)J j are either previously determined 

coolant bulk temperatures or initial input values. The value m^ is the cool- 
ant flow rate and (x.) the mean specific heat of the coolant between B 
and D. Then T (x ) is approximated as 


T (x.) 
f c v l' 


(x AX ‘0 

l t _ 

T »( x -) 
L £ i J 

+ 

j 

T (x. + Ax ) 

f v l i2'J 

J. 

Vi 2 / 
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2 
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(42) 


where the subscript c denotes a calculated value compared with a previously 
determined value or the initial input number. Combining the above three 
equations results in 


T (x.) 
ic v r 




T (x. + Ax 
i v l 


i2^J j 


* r ( x i ) (x.) Ax. 

C p^( x i) cos a ( x i) 


(43) 


This is the internally calculated coolant bulk temperature. For real thrust 
chambers composed of tubes or channels, a cooling efficiency rj should be 
applied to the second term on the right side of equation (43) to account for the 
real geometry effect. Thus, 
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. (44) 
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ic v r 


_ M i+fw^L , 1 r < x i> gw' ( x t ) A y 


-VJ «... — " ■■ — ■'■I' — * 
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4 pt v v v v 


However, the cooling efficiency 17 should be equal to one, if the empirical 
relationship in equation (24) is based upon real thrust chamber data and not 
upon a single tube experiment. 


To start the calculation, a coolant flow temperature distribution must 
be given or approximated to obtain T (x.) through iteration by equation (44) . 

i- C X 

for successive iterations can be obtained internally 


The initial value 
from 


T (X.) 
V i' 

L. _i 



(45) 


Successive iterations are made until the desired convergency is obtained, i. e. , 
the computation is completed when the total heat transfer rate through the 
chamber wall on the gas side SUMQDA in TBL-I and that on the coolant side 
[equation (48)]] (represented by SUMQWI in the present computer program) 
become equal. The specific heat transfer rates through the walls on the gas 

side (q - QW in TBL-I) and the coolant side (q ' = QWI in the present com- 
w w 

puter program) at each section are simultaneously equal. Now the coupling 

of the regenerative cooling cycle and TBL-I [1] is completed. 


SEQUENCE OF CALCULATION 


The numbers of the items that follow in this topic correspond to those 
in Figure 3 a. The calculation sequence at station x = Xj progresses as follows: 

1 . As shown in the flow chart of Figure 3 a, the coolant bulk temperature 
pT (xjl and the gas-side wall temperature distribution [t (x)~| must be 

L i ] o L wg J 0 
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input to initiate the computation, where the subscript 0 denotes the first 
approximated value. The gas-side wall temperature jj w g( x )J ^ is used to 

obtain the heat transfer coefficient on the gas side h (x) at each station 

g 

according to the equation below: 



H 


aw 


T (x 
aw v ’ 


(x) - H (x) 
v ' w 

- pT (x )1 . 

w g J J 


(46) 


where 3 = 0 denotes the first overall iteration loop. Each parameter except 
[T w g(x)J . , on the right side of the above equation, is calculated by the 

equations shown in equations (l) through (19) , or is input. The velocity 
U (x) is the only input parameter in equation (46) which remains constant 

during all iteration at each local station. 

2. The wall temperature on the coolant side T and the heat trans- 

fer coefficient of the coolant flow h^ are calculated by small internal iteration 
loop®, because equation (33) is implicit, 



(equation 33) 


^= 0.025 - 


jjO. 8 p0. 4 


tube 



( equation 24) 


Since each parameter in the previous equations is a function of the axial dis- 
tance x , the argument x is dropped for simplicity purposes. The subscript 
j identifies the iteration number with 3 = 0 indicating the first iteration. 
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is obtained from equa- 


3. 

tion (34) 


The new gas- side wall temperature T 

wgc 


T 

wgc 


h T 
g aw 




T 

W i 


(equation 34) 


where the subscript c denotes a calculated value. The h^ in equation (34) 

It 1 .. The 
[ wgj] 

coolant- side wall temperature T in equation (34) has been obtained previ- 


is still based upon the input wall temperature on the gas side 


w„ 


ously. 


4. The specific heat transfer rate is obtained by any one of equations 
(20), (22), or (23) because of their equivalence represented by equation (32). 
Equation (20) is selected here, 


q' = h (T -T \ 
w g ^ aw wgc J 


(equation 20) 


Another specific heat transfer rate based upon the input gas-side wall tempera- 
ture is obtained from 


h g ( T aw “ [ T wg] j) 


The h term in both equations is based on the temperature [T | When 
g [_ wgj ] 

the overall iterations are completed, the following condition must be satisfied: 


• • j- — 

q = q ’ , because T = T 
w Av Wgc Wg ] 
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5. The coolant bulk temperature has to be corrected at this point by 
considering the heat transferred at section i with the respective x = x^ and 

the input coolant bulk temperatures. 


nn — — “ = — * = m ■■■ ■ ■ - ■ — .... i J ^ yi i ■ — ■■■ 

i c 2 m. C cos o; 

i p£ 

( equation 44) 


Derivation of this equation was shown previously. 

6. New temperature approximations for the bulk coolant and the gas- 
si de wall are predicted, for use in the succeeding overall iterations, from 



and 



T (x) 
wgc 


T (x) 


Xi 


The above procedure from steps 1 through 6 is repeated at each local 
station x = x., and the two total heat transfer rates through the wall are com- 
pared at the end of every overall iteration loop station x = (Fig. 3b). 

A solution is obtained when the two values fall within a small tolerance. 


N ^ N 

2A- .2 A 

i=l i=l 


N 



— Tolerance 
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N 


The expression ^ Q w will be described in equations (47) through (5l) and 


i = 1 


N 


is identified as SUMQWI in the computer program, whereas ^ Q w is based 


i = 1 


upon q w and denoted as SUMQGA. As long as convergence is not attained, 
iterations must be continued with new estimates of 

[V x> ] j ♦ i- 


T„(x) • , 

r 'J + 1 


and 


TOTAL HEAT TRANSFER RATE 

The heat transfer rate through section i , between B and D in Fig- 
ure 4 is Q^x^) according to equation (39). The surface area of this wall 
section is 


2it r(x.) Ax. 

l i 


cos a(x.) 

l 


(47) 


Summation of the heat transferred up to section i = N is equal to 


Z <V • < 48 > 

i = 1 


This amount is the heat which is transferred through the chamber walls 


between x = x, and x = x^„ 

N 

up to X = x^. 


dX N2 

2 


into the coolant per unit time, and not 
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The heat transfer rates through the initial and the final section of the 
contour are those through the area between C and D, and B and C, respec- 
tively, Figure 4. Heat transfer rate through the initial section i = 1 

between x = Xj and x = x t + is 

La 


Q w ( x i) = 


7T r(x!) qJUi) ■ 
cos a (xj) 


(49) 


whereas for the final section at x = Xjxtab 



i '( x ixtab) Sv (^xtab) Ax IXTAB 1 
^“(^xtab) 


(50) 


Integrating the heat transferred from point Xj to per unit time results 

in 


V 


IXTAB - 1 

V X < ) + I VV^^XTAb) 


(51) 


The coefficient r) is used to account for surface area geometry effects; i.e. , 

7] = 0. 5 for double pass cooling if only one path is considered. The above 
amount at x = Xj^tAB mus * : co ^ nc ^ e with the total heat transfer rate calculated 

from the gas phase final iteration value. The latter amount has been denoted as 
"SUMQGA" = M SUMQDA*T 7 n whereas the former, equation (51), will be desig- 
nated as "SUMQWI". These two values are not the same at an intermediate 
section x = x. because "SUMQGA" in TBL-I is the amount between x = x< 

l 

and x=x^, whereas ’'SUMQWI" is obtained between x = x t and 
Ax.„ 

l2 

x = x + . Iterations are performed until the following convergence is 

i 2 

obtained: 
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SUMQGA - SUMQWI 
SUMQWI 


< Tolerance 


SAME DIRECTION COOLANT FLOW 


We have considered the case of the coolant flowing in an opposite direc- 
tion to the combustion products. Now the coolant bulk temperature calculations 
are described for the coolant flowing in the same direction as the combustion 
products. Since rocket nozzles have been built with coolant flow passages in 
either direction and combinations thereof, the up and downstream coolant flow 
simulation of this new concept provides a capability for sectional treatment of 
changing the cooling cycle patterns. Equations in the section entitled Internally 
Calculated Coolant Bulk Temperature which must be replaced for the downpath 
simulation, are shown as follows: In changing the arrow of the coolant flow to 
point in the same direction as the combustion products in Figure 4, the tempera- 
ture of the coolant leaving section i can be determined by 



„ ^ / _ 


m i c p<( x i) 


, (52) 


where the argument x. + — r represents the coordinate at the coolant outlet 

^ 1 u 

location in section i. This equation must replace equation (41). 


The coolant bulk temperature to be calculated at x - x^ is obtained in 


a way similar to equations (42) and (44), 


T ( + *) + - 

fa ( x i - oV 

[ T i ( x i)] i 

T . . i V i 2 ) + ‘ 

2 


T < ( X i) = 

2 


fa ( x , - 1)] i + 

Wi.. 

^w( X i) 


(53) 


l p£ 


( x 0 


. (54) 
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with rj as the cooling efficiency due to geometry effects. Since the coolant 
flow temperature in this case increases toward the nozzle exit, the tempera- 
ture input tables must be arranged correspondingly. 


DOUBLE PASS COOLING 

In carrying out the calculation for a double pass cooling jacket with 
coolant flowing downstream initially and upstream afterwards, we assume, at 
first, that the nozzle wall consists only of down-pass tubes engaged in the 
heat transfer process. A correction is made to the analysis by a cooling 
coefficient r] which represents the surface area exposed to the hot gas covered 
by the downstream cooling tubes, compared to the total surface area. Then, 
the upstream pass calculation is executed in the same fashion neglecting the 
downstream coolant flow part. With each heat transfer calculation process, 
a wall temperature profile is provided. In order to determine the real tempera- 
ture profile for the nozzle wall on the hot gas side, an average from the two 
temperature profiles can be determined. 

The cooling coefficient 17 is usually less than unity for the double 
pass cooling jacket. For coolant flowing in one direction, the cooling coeffi- 
cient may exceed a value of one, since the wall surface area per unit length 
may be greater than the circumferential area due to the ripples formed between 
adjacent cooling tubes. 

In the computer program an option indicator will identify which type 
of coolant flow direction should be considered in the analysis: 

IDUMP = 0 Coolant flow upstream 

IDUMP = 1 Coolant flow downstream 
Modifications made to the existing TBL program are shown in Appendix B. 


EXAMPLE 


In this section of the paper the previously described new concept is 
applied to a thrust chamber nozzle similar to the Space Shuttle's main engine. 
A common chamber down to an area ratio of e = 7 is coupled with different 
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nozzle extensions expanding the combustion products to an area ratio of € = 35 
or e=150 depending on low altitude or vacuum operating conditions, Figure 5. 
The nozzle contours were optimized according to Rao’s method [5,6] to pro- 
vide maximum performance. Since a common chamber, Figure 6, was con- 
sidered for both engines, the orbiter contour had to be modified as indicated 
by the dotted line in Figure 5. In the thrust chamber liquid hydrogen and oxy- 
gen react at a mixture ratio of 6. 0 at a pressure of 3020 psia (212. 33 kgf/ cm 2 ) , 
resulting in a stagnation temperature of 6600°R (3667°K). The free stream 
inviscid flow parameters serving as boundary layer edge conditions such as 

Mach number M , static pressure P , static temperature T and mean 
co co 00 

molecular weight 971, were obtained from the Two-Dimensional Kinetics (TDK) 

computer program [7]. 

First, only the combustion chamber expanding the reaction products 
to an area ratio of e = 7 is considered. In this section the chamber wall is 
regeneratively cooled with liquid hydrogen which flows in an opposite direction 
to the combustion products. The input data for the modified computer program 
are shown in Table 2 and Figure 7. The cross-sectional area variation of an 
individual cooling tube, assumed values for the gas-side wall temperature, 
and coolant bulk temperature as functions of the axial nozzle length, are pre- 
sented in Table 2 and Figures 8 and 9. When a study is performed to optimize 
the cooling jacket geometry, the cross-sectional area in Table 2 and Figure 8 
must be changed in each separate analysis. From such a parametric analysis, 
the best cooling tube geometry can then be selected. In the present example, 
however, the jacket geometry is fixed. Table 3 represents the relationship 
between the specific heat at constant pressure and temperature of the combus- 
tion products in the boundary layer. In order to determine the coolant flow 
heat transfer coefficient, the specific heat, thermal conductivity and viscosity 
for an expected pressure range between 4500 psia and 6000 psia (316. 38 kgf/cm 2 
and 421.84 kgf/cm 2 ) for the coolant fluid must be established as functions of 
temperature. The input data based upon References 4 and 8 are specified in 
Table 4. Additionally required input data can be found in Table 5. The calcu- 
lated temperature distributions on the hot gas-side, liquid coolant-side and the 
coolant are plotted in Figure 10. The total heat transferred through the cham- 
ber wall without considering an enhancement factor is 10 580 kcal/ sec 
(42 000 Btu/sec), whereas the local specific heat flux is exhibited in Figure 11. 
The velocity and temperature boundary layer thicknesses are presented in 
Figure 12 and the momentum and energy thicknesses are plotted in Figure 13. 

The most important result from a performance aspect is the boundary 
layer displacement thickness 6* , Figure 14. This parameter, significantly 
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affected by the wall temperature, reveals by how much the wall contour, iden- 
tical to the inviscid flow border streamline, must be displaced to allow the 
same mass flow condition. A negative sign of 6* means a displacement of 
the inviscid-flow contour towards the thrust chamber centerline. 

If the density across the boundary layer is constant, the profile of the 
mass flow density pu is in principle similar to the velocity profile. Figure 
15a. However, if the density varies the mass flow density overshoots its 
free stream value , especially when the wall is highly cooled, Figure 

15b. The dotted line in either schematic denotes the mass flow density pro- 
file for inviscid flow. Results from the present analysis indicate that the dis- 
placement thickness 6* is negative for the most part of the combustion cham- 
ber to compensate for the strong cooling effect, Figure 14. The performance 
deficiency represented by a thrust loss, Figure 16, down to an expansion ratio 
of e = 7 is already quite large according to the equation [ 1, 2, 3] 


fahrr p' U 1 2 6 cos a 

n 

J exit 

6* 

p ' 

CO 

T r a~ rr 

r oo oo 

e 




^00 CD 


exit 


The corresponding loss in specific impulse is shown in Figure 17. 


To investigate the effect of variable and constant properties necessary 
to calculate the coolant flow heat transfer coefficient, an additional analysis 
was performed using constant values for the specific heat C = 3. 75 Btu/ 

j-lc 

lbm°R, thermal conductivity \ = 0.0000288 Btu/ft s °R and the dynamic 
viscosity = 0. 0000065 lbm/ft s which represents mean values between 


the temperatures of 50°R and 550°R. In comparing the results in Figure 18 
with the ones obtained for variable properties in Figure 10, it is evident that 
the wall temperatures are higher at the throat and lower at an expansion ratio 
of 6=7. This study clearly outlines that most accurate input data must be 
used to perform a reliable analysis. 


Only the chamber section down to an area ratio of 6=1 has been dis- 
cussed. Now, the nozzle extension for the booster engine ranging from an 
area ratio e = 7 to 6 = 35 is treated. For convenience, this nozzle contour 
has been selected, although an analysis for the orbiter nozzle contour would 
be similar. The booster nozzle wall is also cooled by the hydrogen in a 
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double pass cycle. The coolant enters 564 tubes of an area ratio of e = 7, 
flows toward the nozzle exit area (e = 35) and is then turned upstream. The 
wall thickness of each tube varies from 0. 18 to 0. 25 mm toward the nozzle 
exit. All required input data for the downstream and upstream analysis are 
shown in Tables 6, 7, and 8. The resulting wall temperature distributions 
presented in Figure 19 are considerably different for both cooling paths and 
exhibit a minimum in the down-pass section, where the coolant bulk tempera- 
ture reaches a value of approximately 140°K (250°R). At this state the 
hydrogen possesses a maximum specific heat or highest cooling capacity. 

In the real nozzle the temperature differences between the down and up-pass 
cooling tube will come to an equilibrium temperature through lateral heat 
transfer at each local station. Therefore, an arithmetic mean of the different 
temperatures will represent the real nozzle temperature more realistically, 
Figure 20. The individual displacement and momentum thicknesses are pre- 
sented in Table 9, whereas their averaged values are plotted in Figures 21 
and 22. The total performance degradation, expressed in thrust and specific 
impulse loss at the nozzle exit, resulted in AF = 4.742 tons ( 10 470 lbf) 

Jd. Li* 

and A ISP = 7. 687 s (Fig. 23). Heat absorbed by the coolant fluid between 
the injector face and the nozzle exit (e = 35) amounts to 27 000 kcal/s 
(107 000 Btu/s). This method was also applied to identify the area of ice 
formation (wall temperatures less than 460°R) inside the J-2 engine; since 
deposition of ice crystals along the nozzle exit periphery were observed during 
altitude simulation test firings. 1 

CONCLUSION 


A new method has been presented by which the hot gas-side and the 
coolant flow-side wall temperature distributions, as well as the coolant fluid 
temperature variation of a regeneratively cooled thrust chamber, can be deter- 
mined. The analytical formulation is based upon a coupling of the boundary 
layer equations with the heat transfer process through the nozzle wall and 
the coolant flow heat absorption. The new concept has been incorporated into 
the existing JANNAF Turbulent Boundary Layer (TBL) computer program. 

A sample case showing the application of the new calculation process for a 
thrust chamber similar to the Space Shuttle booster engine, has also been 
outlined. Since several empirical relationships such as the friction coefficient 
of the hot gas-side wall, the Stanton number, and Colburn’s equation for the 


1. Analytical Prediction of Ice Formation Inside the J-2 Engine Nozzle 
Contour (200 K Thrust Level) . Memorandum S&E-ASTN-PP (72M-5) 
NASA, Marshall Space Flight Center, January 1972. 
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coolant flow heat transfer coefficient were used and no adjustments for the 
coolant flow turbulence and channel curvature were made, the results are 
only approximate. In addition, this new model could serve as a convenient 
tool for the design of an optimum cooling path and channel geometry concept. 
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Figure 23. Thrust loss due to viscous boundary layer effects. 



TABLE 1. REGENERATIVE COOLING EQUATIONS 


q ' = h (T - T ) 
w g aw wg 


(equation 20) 


h = p U C T — 

g co CO H T 


H - H 

aw w 


'T "I . 

WgJ] 


hi= 0 . 025 _J^ ,.» pM J 

tube Z i \ w. , 


Mr 


(equation 21) 


(equation 24) 


¥ + h if 1 + 7h-^ 

t i \ th 


(equation 33) 


h T + — T 
g aw t w 


h + — 

g t 


(equation 34) 
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TABLE 2. INPUT DATA FOR COMBUSTOR 


1 

X 

Y 

Mach Number 

Pressure 

Static 

Temperature 

1 

-0.950725 

0.802083 

0.220976 

422396.766 

6476.894 

2 

-0.914726 

0.801277 

0.221523 

422122.793 

6476.810 

3 

-0.842728 

0.794779 

0.225972 

421622.242 

6475.952 

4 

-0.770730 

0.781616 

0.235435 

420527.219 

6474.074 

5 

-0.734732 

0.772440 

0.242420 

419691.812 

6472.632 

6 

-0.662734 

0.748573 

0.262271 

417197.340 

6468,317 

7 

-0.590736 

0.716937 

0.292080 

413116.426 

6461.215 

8 

-0. 554738 

0.699832 

0.309721 

410519.758 

6456.660 

9 

-0.518739 

0.682727 

0.329208 

407496.906 

6451.318 

10 

-0.446740 

0.648511 

0.375117 

399765.176 

6437.497 

11 

-0.410741 

0.631406 

0.402511 

394762.754 

6428.419 

12 

-0.374742 

0.614301 

0.433805 

388714.000 

6417.295 

13 

-0.302744 

0.580085 

0.513023 

371954.164 

6385.611 

14 

-0.266745 

0.562980 

0.565406 

359875.809 

6361.932 

15 

-0.230746 

0.545875 

0.628830 

344318.844 

6330.417 

16 

-0. 194747 

0.529072 

0.696138 

327004.102 

6293.651 

17 

-0.158749 

0.515094 

0.766393 

308263.379 

6251.739 

18 

-0.122750 

0.504252 

0.839565 

288302.387 

6204.368 

19 

-0.086751 

0.496347 

0.915600 

267383.355 

6151.280 

20 

-0.057200 

0.491945 

0.980160 

249697.223 

6103.208 

21 

0.000000 

0.488583 

1.193460 

193860.369 

5926.078 

22 

0.003573 

0.488618 

1.227370 

185645.486 

5888.241 

23 

0.007256 

0.488720 

1.260930 

177689.357 

5850.416 

24 

0.011012 

0.488901 

1.293320 

169966.758 

5814.319 

25 

0.014858 

0.489160 

1.325760 

162435.506 

5777.151 

26 

0.018786 

0.489507 

1.358830 

155083.424 

5738.689 

27 

0.022785 

0.489942 

1.392170 

147900.947 

5699.396 

28 

0.030987 

0.491104 

1.459980 

134013.490 

5617.924 

29 

0.035185 

0.491842 

1.494590 

127297.204 

5575.459 

30 

0.043767 

0.493650 

1.565680 

114299.943 

5486.489 

31 

0.048149 

0.494735 

1.602360 

108015.929 

5439.801 

32 

0.057091 

0.497290 

1.678320 

95879.731 

5431.590 

33 

0.061650 

0.498775 

1.717660 

90029.727 

5289.898 

34 

0.070945 

0.502210 

1.799740 

78774.162 

5180.797 

35 

0.080477 

0.506314 

1.887080 

68127.430 

5063.229 

36 

0.095255 

0.513926 

2.005430 

55677.686 

4903.150 

37 

0.110036 

0.522159 

2.016620 

54371.307 

4891.896 

38 

0.125698 

0.531066 

2.035810 

52354.768 

4868.906 

39 

0.150196 

0.545269 

2.056660 

50226.030 

4845.773 

40 

0.186436 

0.566752 

2.087750 

47251.451 

4810.424 

41 

0.222098 

0.588215 

2.108640 

45296.231 

4788.942 

42 

0.257291 

0.609630 

2.133830 

43109.654 

4760.864 

43 

0.292175 

0.631064 

2.157400 

41165.610 

4734.703 

44 

0.326872 

0.652552 

2.183230 

39155.333 

4704.890 

45 

0.361499 

0.674113 

2.206210 

37436.992 

4678.748 

46 

0.396154 

0.695748 

2.228720 

35824.066 

4652.945 

47 

0.430941 

0.717480 

2.251510 

34260.629 

4626.427 

48 

0.465973 

0.739363 

2.274500 

32744.594 

4599.258 

49 

0.501365 

0.761452 

2.297430 

31295.269 

4571.973 

50 

0.537202 

0.783776 

2.321400 

29851.381 

4543.057 
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TABLE 2. (Continued) 


1 

X 

Y 

Velocity 

Molecular 

Weight 

Coolant Area 

Coolant 

Temperature 

WaU 

Temperature 

1 

-0.950725 

0.802083 

1175.830 

13.590460 

0.000053330 

460.000 

1360.000 

2 

-0.914726 

0.801277 

1178.720 

13.597353 

0.000053330 

450.000 

1360.000 

3 

-0.842728 

0.794779 

1202.300 

13.597646 

0. 000053330 

430.000 

1360.000 

4 

-0.770730 

0.781616 

1252.440 

13. 598282 

0.000053000 

410.000 

1360.000 

5 

-0.734732 

0.772440 

1289.430 

13.598768 

0. 000052500 

398.000 

1360.000 

6 

-0.662734 

0.748573 

1394.480 

13.600157 

0.000051100 

380.000 

1360.000 

7 

-0.590736 

0.716937 

1551.980 

13.602585 

0.000049800 

366.000 

1360.000 

8 

-0.554738 

0.699832 

1645.040 

13.604074 

0.000049100 

350.000 

1360.000 

9 

-0.518739 

0.682727 

1747.700 

13.605855 

0.000048500 

340.000 

1360.000 

10 

-0.446740 

0.648511 

1988.940 

13.610415 

0.000047100 

320.000 

1360.000 

11 

-0.410741 

0.631406 

2132.440 

13.613397 

0.000046500 

312.000 

1360.000 

12 

-0.374742 

0.614301 

2295.920 

13.617002 

0.000045800 

303.000 

1360.000 

13 

-0.302744 

0.580085 

2707.390 

13.627159 

0.000044400 

287.000 

1364.000 

14 

-0.266745 

0.562980 

2977.420 

13.634497 

0.000043700 

280.000 

1368.000 

15 

-0.230746 

0.545875 

3301.920 

13.646259 

0. 000043000 

270.000 

1370.000 

16 

-0.194747 

0.529072 

3643.090 

13.658344 

0. 000042400 

262.000 

1372.000 

17 

-0.158749 

0.515094 

3995.360 

13.672035 

0.000041700 

255.000 

1374.000 

18 

-0. 122750 

0.504252 

4357.770 

13.687345 

0.000041000 

247.000 

1380.000 

19 

-0.086751 

0.496347 

4729.150 

13.704345 

0.000040868 

240.000 

1380.000 

20 

-0.057200 

0.491945 

5040.040 

13.719516 

0.000040868 

232.000 

1375.000 

21 

0.000000 

0.488583 

6035.780 

13.765060 

0.000040868 

222.000 

1320.000 

22 

0.003573 

0.488618 

6186.780 

13.768245 

0.000040868 

221.000 

1310.000 

23 

0.007256 

0.488720 

6334.890 

13.771604 

0. 000040868 

220.000 

1300.000 

24 

0.011012 

0.488901 

6473.390 

13.790009 

0.000040868 

219.500 

1295.000 

25 

0.014858 

0.489160 

6612.360 

13.799884 

0. 000040868 

219.000 

1285.000 

26 

0.018786 

0.489507 

6752.420 

13.809909 

0.000040868 

218.000 

1280.000 

27 

0.022785 

0.489942 

6891.970 

13.820396 

0. 000040868 

217.000 

1272.000 

28 

0.030987 

0.491104 

7170.700 

13.842066 

0.000040868 

216.000 

1255.000 

29 

0.035185 

0.491842 

7310.530 

13.852492 

0.000040868 

215.000 

1250.000 

30 

0.043767 

0.493650 

7592.250 

13.872702 

0.000041000 

214.000 

1235.000 

31 

0.048149 

0.494735 

7734.490 

13.883344 

0.000041100 

213.000 

1225.000 

32 

0.057091 

0.497290 

8022.470 

13.905751 

0.000041400 

211.000 

1220.000 

33 

0.061650 

0.498775 

8168.540 

13.915809 

0.000042000 

210.000 

1210.000 

34 

0.070945 

0.502210 

8456.700 

13.936370 

0.000042400 

209.000 

1200.000 

35 

0.080477 

0.506314 

8770.780 

13.958197 

0.000043200 

207.000 

1188.000 

36 

0.095255 

0.513926 

9166.640 

13.987167 

0. 000044000 

205.000 

1180.000 

37 

0.110036 

0.522159 

9205.220 

13.993980 

0.000045200 

202.000 

1170.000 

38 

0.125698 

0.531066 

9270.930 

13.996067 

0.000046000 

199.000 

1155.000 

39 

0.150196 

0.545269 

9342.800 

14.000566 

0.000047700 

194.000 

1145.000 

40 

0.186436 

0.566752 

9449.400 

14.003944 

0.000050000 

188.000 

1145.000 

41 

0.222098 

0.588215 

9522.320 

14 . 006823 

0.000052500 

182.000 

1160.000 

42 

0.257291 

0.609630 

9607.770 

14.009449 

0.000055000 

178.000 

1180.000 

43 

0.292175 

0.631064 

9686. 870 

14.012815 

0.000057200 

172.000 

1180.000 

44 

0.326872 i 

0.652552 

9771.880 

14.015822 

0.000059400 

167.000 

1170.000 

45 

0.361499 

0.674113 

9847.080 

14.018885 

0.000061800 

163.000 

1164.000 

46 

0.396154 

0.695748 

9919.870 

14.021893 

0.000064000 

160.000 

1150.000 

47 

0.430941 

0.717480 

9992.570 

14.024915 

0.000064000 

155.000 

1140.000 

48 

0.465973 

0.739363 

10065.200 

14.027241 

0.000064000 

152.000 

1125.000 

49 

0.501365 

0.761452 

10136.800 

14.029441 

0.000064000 

150.000 

1113.000 

50 

0.537202 

0.783776 

10210.600 

14.031649 

0.000064000 

144.000 

1100.000 


45 
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(The enhancement factor 77 is assumed 1. 0 ) . 



I 

Specific Heat 
(Btu/lbm s) 

1 

0.6199999973 

2 

0.6550000012 

3 

0.6799999997 

4 

0.6950000003 

5 

0.7049999982 

6 

0.7199999988 

7 

0.7282169983 

8 

0.7282169983 

9 

0.7282169983 

10 

0.7282169983 

11 

0.7282169983 

12 

0.7282169983 

13 

0. 8833189979 

14 

0.8843249977 

15 

0.8854160011 

16 

0.8893399984 

17 

0.8902480006 

18 

0.8906119987 

19 

0.8907269984 

20 

0.8920999989 


Temperature 

(°R) 

400.000 

800.000 

1200.000 

1400.000 

1600.000 

2000.000 

2500.000 

3000.000 

4000.000 

5000.000 

5850.000 
5926.078 
6103.208 
6151.280 
6204.368 
6403.409 
6451.318 
6470.760 
6476. 894 

8000.000 









TABLE 4. PHYSICAL PROPERTIES OF LIQUID HYDROGEN 
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(Regenerative cooling in the opposite direction; Injector to e = 7) 



INPUT DATA OF SSME BOOSTER NOZZLE (DOWN PASS) 
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(SSME Booster Engine Double Pass Cooling from e = 7. 0 to 35. 0) 
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(The enhancement factor is assumed 1.0) 




TABLE 8. INPUT DATA OF SSME BOOSTER NOZZLE (UP AND DOWN PASSES) 
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TABLE 9. CALCULATED DISPLACEMENT AND MOMENTUM THICKNESSES ALONG NOZZLE WALL 
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DESCRIPTION OF PROGRAM INPUT 


MZETA = n 
IPRINT 

IXTAB 
ICTAB 
ITWTAB 
TO = T 0 
P0= P 0 
GAM.0 = yQ 
ZM.U0 = M 0 
ZMVIS 
ZNSTAN 
DXMAX 
THETAI = 0. 

l 

PHII = 0 . 

l 

EPSZ 


Input Data 

Exponent in velocity profile power law 

Print option at every calculated point (= l) 
or at input intervals (= 0) 

Number of points in x = f(y) and x = ) tables 

Number of points in = f(T) table 

Wall temperature option = 1 (must be input) 

Stagnation temperature, °R 

Stagnation pressure, lbf/ft 2 

Stagnation specific heat ratio 

Stagnation viscosity, lbm/ft-s 

Exponent of viscosity temperature law 

Boundary layer interaction exponent 

Maximum step size 

Initial value of momentum thickness, ft 

Initial value of energy thickness, ft 

Geometry option - Axisymmetric = 1. 

Plane = 0. 


RBAR 


FJ = J 


G * g 


Gas constant at stagnation, ft-lbf/°R-lbm 

Conversion factor between thermal and work 
units = 778. 2, ft-lbf/Btu 

Acceleration of gravity = 32. 174, ft-lbm/lbf-s 2 
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SCALE 

ITZTAB 

IDUMP 

FLOWRT 

MASSL = m^ 

RAMDW = X 
w 

COEFCL = 7 ] 
TUBEN 


(i) 

Specific 

(ii) 

XITAB 


YITAB 


ZMTAB 


PE TAB 


TETAB 


UETAB 


SMTAB 


ALTAB 


TLTAB 


Contour scale factor 

Number of points in temperature versus C^ , X^ 
and \i^ table 

Coolant flow option 

Same direction = 1 
Reverse flow = 0 

Combustion chamber mass flow rate, lbm/ s 
Coolant mass flow rate, lbm/ s 

Thermal conductivity of the chamber wall, Btu/ft-s°R 
Cooling coefficient (surface area effect) 

Number of cooling tubes 


Input Tables 


C (CPTAB) versus temperature T (TITAB) 

P 

Axial distance, ft 
Radius, ft 

Mach number at boundary layer edge 

Static pressure P^ at boundary layer edge, lbf/ft 2 

Static temperature T at boundary layer edge, °R 

Velocity at boundary layer edge, ft/s 

Mean molecular weight 311 at boundary layer edge 

Cross-sectional area of each cooling tube, ft 2 


Assumed coolant temperature 


, °R 


60 



TWTAB 

Assumed wall temperature on the gas-side TT ] , °R 

L W SJ o 

THITAB 

Wall thickness of the cooling jacket, ft 

(iii) TZTAB 

Coolant temperature table used to obtain C . , A 
and n , °R ** 

CPLTAB 

Coolant specific heat C . , Btu/lbm-°R 

p£ 

RAMTAB 

Thermal conductivity of coolant , Btu/ft-s°R 

ZMYTAB 

Viscosity of coolant , Ibm/ft-s 

DESCRIPTION OF PROGRAM OUTPUT 

The following parameters are printed out in addition to the original TBL 
computer program results [ 3] : 

RBAR =^/3R 

Specific gas constant, ft-lbf/lbm°R 

PRANDT = Pr 

Prandtl number of the free stream 

GAME = y 

oo 

Specific heat ratio at the boundary layer edge 

SMOL = 9H 

Mean molecular weight, lbm 

COSAL = cos a (x) 

Cosine of the wall angle 

DELFA = AF d t 

B, L, 

Thrust degradation due to turbulent boundary 
layer effects downstream of the throat only, lbf 

THRUST = F 

Vacuum thrust, lbf 

DEFTHR = AF/Fx 100 Percent of thrust degradation 

TBLISP = -A I 

sp 

Specific impulse loss due to turbulent boundary 
layer effects, s 

THRUSA 

Thrust at sea level, lbf 

VISP = I 

sp 

vacuum 

Vacuum specific impulse downstream of the 
throat only, s 
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AISP = I 

Sp sea level 
DMASSL = p 

HL = h 

i 

QWI = 

REYL = R 

e i 

SUMQGA 

SUMQWI 

TEMPRL = T /T. 

W i 1 


Specific impulse at sea level of the throat only, s 

Mass flow density of the coolant fluid, lbm/ft 2 -s 

Heat transfer coefficient of the coolant fluid, 

Btu/ ft 2 - s° R 

Specific heat transfer rate based upon calculations 
for the coolant flow side, Btu/ft 2 -s 

Reynolds number of the coolant fluid based upon 
tube diameters 

Total heat transfer rate, Btu/s 

Total heat transfer rate (includes cooling flow 
calculation), Btu/s 

Temperature ratio 


TLCA = T^ Calculated coolant temperature, °R 

TWGCA = T Calculated wall temperature on the gas side, °R 

wgc 

TWL = T Calculated wall temperature on the coolant side, °R 

W .f 

DIATUB= 2 \l Af.ykg/ 71 ’ Equivalent diameter of the cooling jacket, ft 

THICK = t Chamber wall thickness (input value), ft 

George C. Marshall Space Flight Center 

National Aeronautics and Space Administration 

Marshall Space Flight Center, Alabama, February 11, 1972 
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APPENDIX A 


DERIVATION OF EQUATIONS (33) AND (34) 


Using equations (20) and (22) with equation (32), we obtain 


g 


(T - T ) = - 7 — (T - T ) 
v aw wg 7 t wg w, 


wg 


Rewrite the above equation, as 


m 

, m W T 

h T + — — w 

g aw t I 


wg 


h + 
g t 


w 


(equation 34) 


Equations (20) and (23) reduce to 


h (T - T ) = h /T -T\ 

g' aw wg 7 i I w^ i 1 


so that 


T = T ~ T c\ 

wg aw h g \ ij 


Substitute equation (34) into the above equation, then 


w 

h T +-— T 
g aw t w 


h + 

V g t 


w 


(V T ') 


t It - t 

aw h g \ ^ i 
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APPENDIX B 

TBL MODIFIED COMPUTER PROGRAM LISTING (TBLREG) 
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TPFS. BARCON 


BARCOOO 1 


C 

c -- 

c 


c 

c 

c 

c 

c 

c 

c 

c 

10 

1 

1 1 


SUBROUTINE BARCON 


BARCON -- CONTROLLING SUBROUTINE B ARCQ002 

COMMON /COOL/ I C OOL . I D UMP . I T Z T AB . AL . COEF CL i CPL . OEL XB A , D 1 A TUB . /COOL/ 

1 FLOART,MASSL.PRANOL.RAMDL,RAMO*,REYL,SUMQGA.SUMe*l , /COOL/ 

2 THICK, TLO, TLl .TL2.TLCA.T0LITE,TUBEN,TWGCA,ZMYUL. /COOL/ 

3 CPLTAB120) .RAMTAB120) .TZTAB120) ,ZMYTAB(20l , /COOL/ 

H ALTABI ICO) ,THITAB( IOC) .TLCTABUOO) .TLTABI 1001 . /COOL/ 

& Ti&GTABI 1001 /COOL/ 

REAL MASSL /COOL/ 

COMMON /INPUT/ I OX M A X , I C T AB , 1 PR I N T , I T* T AB • 1 * T AB . MZ E T A . 0 X M A X , /INPUT/ 

A EPSZ.FJ.&.GAMO.PO.PHI I .PIE.PRANOT.RBAK.SCALE.TO. /INPUT/ 

K THETAI .TOLCFA.TOLZET.TOLZME.ZMUO.ZMVIS.ZNSTAN /INPUT/ 

COMMON /INTER/ CFA6T .cPa6P*CHPAR 1 .D* •DXRHO.HE.Hi*. I bEG.MZETAM, /INTER/ 

A O0MZeT,Ph1P.PREI03,RH0E»RH0UE,RMZETA,ThETAP, /INTER/ 

K xIBASE.XIENO.zETaTM.ZMzETA.ZMzETm.ZMZETP /INTER/ 


COMMON /LOOKUP/ 

1 

2 


1CX.1MX,1PX,[RX,ISX.ITP0S,IT«X.ITX,1UX.IXP0S.IYX, 
1 ZX ,CC X ( 6 ) , CMX ( 6 ) , CPX I 6 ) , CRX ( 6 > .CSX 1 6 ) «C TWX ( 6 > . 
CTXI6) . C UX ( 6 1 , C Y X ( 6 ) » C Z X ( 6 ) 


/LOOKUP/ 

/LOOKUP/ 

/LOOKUP/ 


COMMON /NHANCE/ I E X , CE X ( 6 ) , E NH T a B ( 1 00 ) 


/NHANCE/ 


COMMON /OUTPUT/ 
A 
K 


BOELTa ,CF ,CH .DELTA .DELSOT .OELSTR ,FLaT .FORCE , HG , 

PE,PH1.U».SUM80A,TE.THETA.T*,UE.X.XLARC,YR.ZI. 

ZZ.Z3.ZN.ZB.ZETA.ZME 


/OUTPUT/ 

/OUTPUT/ 

/OUTPUT/ 


COMMON /TABLES/ PE T AB ( 1 00 ) , SM T AB ( 1 00 ) , TE T AB ( l 00 ) , T* T AB ( 1 00 I . /TABLES/ 

1 UETAB1 100) .XITABI 100) .Y1TABI 100) .ZMTABl 100) /TABLES/ 


U I mens I ON DPHIRK(r) .OTHERKIM) ,XCCP( 100 I .YCCPI 100) 

IF 1 1 COOL »EU. 0) GO TO H 
ITER • o 
ITER ■ ITER ♦ 1 
ARITE U,l> ITER 

format t ihi .ssx.&nhrEgEnerative cooling *all temperature iteration 
I NUMBER , I 3////////// ) 


MZfcTAM « MZETA - l 


ZMZETA-MZETA 

BARC003I 

ZMZETP-ZMZETA+l . 

BARC0032 

zmzetm*zmzet a- 1 • 

0AKCOO33 

RMZETA- l ./ZMZETP 

BARCQQ3H 

OOMZET- l ,/ZMZETA 

BARC0035 

X»XITAB( 1 ) 

BARC0036 

DX*0 • 

B A RC 0037 

XL A KC*0 » 

B ARC0038 

sumqua-o. 

BARC003? 

SUMw^A ■ 0*0 


SUM W rt I * 0*0 


FQRCE-O. 

B ARCOOHO 

FLAT«0. 

BARC00H1 
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r> r> n 


BARC00N2 
B ARC00H3 
B ARCOONN 


QW 

hg 

I XP05* 1 
ITPOS « J 
ICX » 0 
IEX - 0 
IRX a 0 
1SX « 0 
IUX a 0 
IZX - o 
IHX*0 
ITX«0 
IPX-0 
IYX-0 


0.0 

0.0 


IT#X»0 

OxRHO-O. 

IBEe a 2 

CFA6T a ,002 

l ST ART a 0 

IF (THETAI .LE. 0.0) 60 TO 2 

ZETA a <PHII/THETAI) ♦•RMZET A 

60 TO 3 

2 CALL START 

ISTART a | 

3 CFA6P a C F A 6 T 

IF (ICOOL .EO. 0) 60 TO << 


OELXQL a 0.0 - 

DELXNE a ABSIXITABI2) - XITABllI) 

OELXBA a (OELXOL * DELXNE) /2.0 
AL a ALTAB(I) 

TLI ■ TLTAB ( 1 ) 

THICK a THITAB(I) 

TL2 a TLTABI2) 

TLO a TLI 

CALL XNTERP (TLI »ZHYUL»ZP. IZX. TZ TAB t ZMY T AB , I TZT aB , C ZX , ITPOS) 
ITPOS a )ZX 


PIATUB a 2.0*S«BT(AL/PIE) 

REYl a HASSL*DIATUB/(AL»TUBEN»ZMYUL) 

H PHI a PHI I 

THETA a THETAI 

x i base a xitabii) 

XiEnD a X I TAB ( ) XT AB ) 

IF UXTaB .LE. 1) AO 10 IS 
DXRHO a ( X I T AB ( 2 ) “ XlBASE)/l0.0 

IS CALL BARPRO(l) - 

CALL BARPRO(S) 

TAG T AB ( l ) ■ TKGCA 
TLCTAB(l) • TLC A 

CALL XNTERP ( X, YR, YRP. tYX. X I TAB . Y l T AB , IXTAB, CYX, IHX ) 


save initial y and delstr. 

DEL a DELSTR 
YMlN a YR 
C 

ONOC • SORT ( I. ♦ YRP • YRP ) 


B ARCOOHS 
BAKC00N6 
BARC00H7 
BARCOONB 
BARC00N9 
BARCOOSO 
BARC00S1 
BARC00S2 
BARC00S3 


BARC00S6 

BARC00S7 

BARC00S8 

BARC0059 


BARC0062 

BARC0063 


BARC0066 

BARC0067 


BARC006B 

BARC0069 

BARC0070 

BARC0071 

BARC0072 

BARC0073 

BARC007M 

BARC007S 
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XCCPU) ■ X ♦ 0EUST8 • YRP / ONOC 

YCCP(l) a YR - OELSTR / ONOC 

IF uxtab *le. n return 

00 20 I a IBEG.IXTAB 
XNE«aXlTAB( I I 

IF (ICOOL • E6 • 0) GO To 16 
AL a ALTABI I I 

thick * thitabi i i 

OELXOL » ABS ( X I TAB ( 1 ) - XITAB(I-ll) 

IF (I ,GE. IXTABI GO To 3000 
OELXNE a ABSUITABUMJ - AITABtll) 

TLO a TLTABI I-n 
TL 1 a TLTAB ( I ) 

TL 2 » TL TAB ( 1 * I ) 

GO TO 3001 

3000 OELXNE ■ 0*0 
TLO » TLTAB (1-1) 

TL1 ■ TLTAB ( I ) 

TL2 a TLl 

3001 OELXBA a (OELXOL ♦ 0ELXNE)/2.0 
16 X M A 6 a (ABS(XnEw) ♦ ABS1XII/2.0 

OX iNTaXNEW-X 

NX ■ OXINT / 0 XM A X + 0.9? 

IF (NX .GT. 0) GO TO 18 
NX » 1 
18 ZNXaRX 

OXaOXJNT/ZNX 
0 X02*0 X / 2 . 

OXRHOaDX/lO, 

DO 30 1 N X ■ I .NX 
PHlOLOapHl 
THEOLDaTHETA 
XOLOaX 

OPH-lNKt l laDXapHIP 

0THERK1 1 )aOX»THETAP 

XaXOLD*OX02 

DO NO l RK*2 » N 

IF (IRK .NE, N) GO TO N N 

X a XOLD * DX 

IF ( ABS ( ( X - XNEkI/XMaG) .GT. 1.0E-6) GO TO H3 
X • XNEw 

N 3 PHI « PHIOLD ♦ OPHlRK(lRK - 1) 

ThETAbTHE0LD*DTHEKK(1RK-1) 

GO TO <45 

N N PHI ■ PHIOLD ♦ OPHIRKIIRK - llaO.SO 
ThETA»ThE0L0+DTHERK(IRK-1)*.S 
NS IF (PHI ,LE. 0.0) GO TO 62 

if itheta .le. o.o) go to 62 
call barproc irk i 

DPH IRK ( I RK ) aOX»PH I P 
NO DTHERK ( 1 RK I ■ DX»THETaP 

PHlaPHI0L0*(DPHlRK(I)*2.*0PHlRK(2I*2.*0PHlRK(3)a0PHIRK(Nl)/6. 
THETAaTHE0LD+(DTHERKllU2.»DTHERK(2)*2.*DTHERK(3>»DTHERK(N))/6. 
IF (PHI .LE. 0«0) GO TO 62 
IF (THETA .GT. 0.0) GO TO 72 
62 *RITE(6,63) X, ZME, THETA, PHI 


BARC0076 

BARC0077 


BARC0080 


BARC0082 

BARC0083 


BAKC0086 

BARC0087 

BARC0088 

BARC0089 

BARC0090 

BARC0091 

BARC0092 

BARC0093 

BARC009N 

BARC009S 

BARC0096 

8ARC0097 


BARC010N 
BANCO 1 05 

BARC0108 


BARC01 13 
BARC011N 

BARC01 17 
BARC01 18 


BARC0I21 
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n n n rs n r\ n n n n o r> r> o n 


H1H **BARCON F A I LURE * * • AXIAL DISTANCE X * . I PE M * 7 • BARC0122 

SX, 11HHACH NO* * * EM*7t 2X» 6HTHETAI «. E i 4 • 7 . 2X» BARC0123 

6HPHII *. E 1 H • 7 / H4H THETA I OR PHU COMPUTED AS NEGAT l VEB ARCO 1 2M 
/ 64H *CHECK CONTOUR AND MACH NUMBER DISTRIBUTION T A BL E SB ARCO 1 25 
H FOR ERRORS./ HQH *M0RE INPUT POINTS MAY BE REQUIRED TO ADEQUATEBARCQ126 
5LY OESCRIBE DERIVATIVE VALUES ALONG THE CONTOUR AT THIS POINT. / BARC0127 
6 * A SMALLER RUNGE^KUTTa STEP S I z£ MAY BE REQUIRED TO AMF QH A Tfi A PC 0 l 28 


63 FORMAT < 
1 
2 

3 OR ZERO 


J fTFCp/T I p p VAl t*Ff 


rri y appfox tmate 

wi»U. to m ? M C , L i 
CALL QUITS 
72 CALL BARPRO ( 1 ) 


SELECT minimum Y and its corresponding oelstk. 



IF < Y 

N.GT.YM1 

N) 

GO TO 

29 




DEL 

* DELSTR 






Y M I N 

S YR 






29 

IF < 

IPRINT 1 

>LE« 

* G) 

o 

►— 

o 

30 



CALL 

BARPRO { 

5) 





30 

CONT 

INUE 







call 

XNTERP 

( ) 

U YR. 

YRP . 

IYX 

i XI 


onoc 


m 

SwRT t 

1 • ♦ 

YRP 

• Y 


xccp 

( I ) 

m 

X ♦ DELSTR 

♦ Y 

RP / 


YCCP 

( I ) 

a 

YK - 

UELSTF 

i / 

ONOC 


IF ( 

IPRINT * 

i G T i 

> 0) 

GO TO 

20 



call 

BARPRO 1 

IS) 






20 
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yitab* ixtab* cyx, I M X ) 


ip ) 
ON 0 C 


T WGT AB t I ) 
TLCTAB t I I 
CONTINUE 


TWfiCA 

tlca 


Ym In 

DEL 

RPOT 

RPOT 


HJNIMUM Y VALUE FOR NOZZLE. 

DELSTR CORRESPONDING TO MINIMUM Y 
ThE POTENTIAL ThROAT RADIUS* 


t T hROaT ) • 


YMIN 


- DEL 


*RlTE(6,l COO I RPOT 


normalize Table uf corrected 
throat radius* 


contour points using the potential 


xccpi i ) 

YCCpt i ) 
DO 7V I 
XCCP( I ) 
YCCpt I > 


XCCPU) / RPOT 
YCCPU) / RPOT 
I BEG • i AT AB 
XCCP(I) / RPOT 
YCCP ( I ) /RPOT 


BS 
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*R I TE ( 6 t 1 CO 1 ) 

IF IISTaRT *LE • 0) GO TO 65 

WRITE C A. I 010) XCCPI 1) » YCCPt l ) • ( I .xCCpt I ) .YCCP ( I ) . I * 1 BE 6 . I XT AB ) 
GO TO 86 

ARITE<6,102C> ( 1. XCCPI II • YCCPII). I » It IXTAB ) 

IF (ICOOL *EQ. C) RETURN 

IF ( ABS( (SUmqdA*COEFCL - SUMQW I ) /SUMQWI ) .LT. TqLITE) RETURN 


mrcouv 

B A wCu i JU 
BARC0I31 

BARC0134 
B ARCO I 35 
B ARCO 136 
BARC0137 
BARC01 38 
8ARC0139 
B ARCO I NO 
8 ARCO 1 H2 


BARCOMS 
BARCOl H6 
BARC01N7 
BARCO MB 
8ARC0M9 


BAKC0153 
BANCO l SH 
BANCO 1 55 
BARC0157 
BARCO 1 56 
BARCOl 59 
B ARCO 1 60 
BARCOl 6 1 
BARC 0162 
BARCOl 6 3 
BARC016H 
BARCOl 65 
BARCO 166 
BANCO l 67 
BARCOl 66 
BARCO 169 
BARCOl 70 

BARCO 1 73 
BARCOl 7 4 


BARCO 179 
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DO 87 J_ - i,lXT_A8 

TRTAB(I) ■ (TWTAB(I) ♦ TWGTAB ( I ) ) /2.0 
87 ILIABCU -•-lJLTA8tI> ♦ TLCTAB II >» / 2 *0 
GO TO 10 

1000 FORHaT-H HI » 2?X IHThRO AT RADiUS CORRECTED FOR DISPLACEMENT , BARCOlBl 

1 UhTHICKNESS ••IPE1S.8//1 BARC0182 

1001 F0RMATlLM0*29X.H8HTA8tE OF NORMALIZED ConToUR PoInTS CORRECTED F ORBARCO 1 83 

1 23H DISPLACEMENT THICKNESS // 37X.10H0ATA PO I NT , 1 OX , l HX , 2**X , 1 H Y/ » BARCO 1 8H 

1010 FORMAT ( *tQX « BHM ■ l.i MX. 1PE1S.8, 10X, E1S.8 / C BOX. IS. BARC018S 

1 SX, E1S.8. 10X. E1S.8 ) ) BARC018B 

lQzO.filRMAT.t BOX. IS* SX. IPElS.S. lOX. E|S.B ) BARC0I87 

END BARCO 1 88 
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SUBROUTINE BARFKO fTNDr 


TOMliON TcDf 1 IF/ If JNT »AFInT»BF1NTiCF INT *TF iixT /con if/ 

COMMON 7c00L(T"iC0U L , lOUMP. I TZTabYa L » C 0 E FC L • C r L • t) E L X B A iD f A tub • /COOL/ 

1 f LOWRT^ASSLiPRANDloRAMULiRAMDw.RfcYLfSUMyGA iSUMQKl , /COOL/ 

2 TM 1 C k • T L 0 • TL l » T L 2 • T L C A , T G L I TE, I Ufl&N , T *GC A » ZM YUL . /COOL/ 

J CPL TaB l 20) ,KAMTAB (20) • TZTABi 20) , ZMYTAB < 20) • /COOL/ 

H- - AL T Aa I r oori rnrtTA^UOO^ ,TLCTAB( IOC) ,TLTAtiUOQ| • /COOL/ 

b T*GT A B<iOG) /COOL/ 

REAL MASSL /COOL/ 

C U HHO N / C b E V A L / NU C T AU , i S t KO J i T JCTCUG , GHITT2 twTO C Ml ,F0 h A X ,CP OiW, /CSCVTL/ 

A S0,TCTAB(20I ,CPTabUO) i6CP(2o) .CCP120) ,0CP(2O) # /CSEVAL/ 

K 6TAB(2d» iHTABr2C) iBARbi < 20) ,t>ARB2 ( 20 ) .BARB3 I 20) /CsEyAL/ 


Common /input/ Toamax # 1ctab. i print ,it*tab,i at ab ,mxet a ,dxmax • 

A EP5j:iFJi0t6ANOiP0»PHlI • P 1 E^PR ANoT^KbAR • SCALE_^TO^ 

~K THEjAi »TOLcFA, TOLlET » TO LZME , Z«Ui&* AM V I 5 * ZNSTAN 


/INPUT/ 

/INPUT/ 

/INPUT/ 


COMMON /InTEk/ O a g T * c F a G r • c HFXRT » D X * DXKHO f HttHifr»lBEG>HZETAMi /INTER/ 

A OOMiET,PHIP»PRE103>RhOE,RHUUE,RM2ETAiTHETAP, /INTER/ 

~K XTBA^nENOiZETATh.ZMiETA^ZhAETM./MZETP /INTER/ 


C 


c 

L 


C 

c 


c 


TT 


COMMON /LOOKUP/ ICX » I MX . I PA # J KA V ISA • iTFOS » I T«X » i T A i i UX*iXp65> I YX * /LOOKUP/ 


1 i 2 X » CC X ( 6 ) iCMX(6)tCPX(6) • CR X t 6 ) » C b A ( 6 ) tCT*Xl6) 9 /LOOKUP/ 

1 l C I X t 6-TTCuX VC yx u r» CZX ( 6 ) /LOOKUP/ 

COMMON /NHaNCE/ IExiCETUI • ENH T AB ( 100) ” /NHANCt / 

CUM ff ON /OUTPU T / BU E LTA,CF>Cff t bELTA~,FE lSOT , _ D El S T ft » F L A T , F QRCE~tHG, /OUTPUT/ 

A PE, Phi .WW.SUMQOAiTEi THETA,T*,UE* A,XLARC t YR # 2l , /OUTPUT/ 

K Z27ZX, 7R aS . zETA.ZME /OUTPUT/ 

COMMON /SAVEU/ A.B.CiZT l,ZllPfZl2.Zl2PiZI3 t Z|3P#Zi<*,Zib.ZJ6.2I7 /SAVED/ 

— COMMON /TABLES/ PLxAb ( 1U0 J • SMT AB U DO « • I ETaBT lUCT) • T AT A B (TOO ) # ” /TABLE 57” 

I UEyABi 1U0) «XI TABU 00 I i Y I TAB ( |0Q) iZmTAcU I 00) /TABLED/ 

DIMENSION ZlNTPRUo* 


OATA ( Z INTPTT TIT » 1 • l fTCT 75HZT4 *#6HZIS -*68X10 a t*HZl7 »• 
A 6 HZ I j P » * 6HZ I | ■ 1 68Z l 2 «,&H2I3 -* 6HZ i 2 P * 1 6 H Z 1 JF ■ / 


GO TO ( H,M|3,H|S) * I NO 

CALL XNTERPrr. ZME t XH CFrUMlU X i T A B t ZMTAB » iXTAB, CMX. iXPOS) 


I APOS- I MX 

CALL XNTERr (X.TETrEPilTXtXITAB.TETAB.lXtAB.LTX. IXPOS) 

CALL XNTERP «X,PEipEP.lP‘X.XlTAB,PETAtt # iXTAtt,CPX»iXPGS> ... 

CALL XNTERP (X.IiEiuEPilVxTXITAB.UETAB.lXTAB.LUX.lXPOSi 

Call XNTERP (X»SMO L »SMQLP* ISX.XlTABtSMTABt lXlAb»CSX» IXPOS) 

cael^^sevali i *te iCPeime) - 

RBAR ■ 15HS*Q/SM0L 
ROV ■ R B A R / F J 

G A mE ■ CPE/ (CPE - rOjJ 

PRANOT ■ H. 0*UAME/ ( 9* 0«(»AME - S«Q) 

UE202 » UE»UE/2*0 ... r 

HEP»FJG*CPE«TEP 

rhoe.pe/te/kbar 

• of A G HOST l C* ThETEST FOR EG U A E j T Y BETWEEN NON- INTEGERS M A Y NOT BE MEANINGFUL* 

~ I F (QXRHO » NE ♦ 0*0) GO TO 2Qj . 

NHOEP ■ 0*0 

GO To 210 

Zt)T rr T x ♦ GT7 XTbXSZT gTTTO 203 
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21 * RHOE 

"Z 3 * l . 

CO TO 20H 

203 CALL XNTERP ( X-D XR h 0 , ZH , Z HP , 1 P X , X I T A b , PET AB » 1 X T Ab » CP X , 1 XPOS > 

CALL XNTERP (X-0XR H 0.Z5»2SP,lTX,XlTAb,TETAB,lXTAB,CTX,lXPaS» 

Call xnteR'p ( a-uxKh°*smi ,smip, i sx , x i t as ,smt ao , i *t ab , cs a , i xpi>s ) 

R 1 * 15H5.0/Srti 
Zl » ZR/Zs/R 1 
Z 3® • 5 

ZQH IF (X .LT. X I E NO ) CO TO 206 
Z 2 * RHOE 
Z3»l , 

CO TO 207 

206 CALL ANTERP ( X .0 X R m 0 . Z H , Z H P . I P X , X l T Ab , PE T aB i 1 X T AB . CPX . I XPOS ) 

CALL XNTERP iX^0XR H 0,Z5,Z5P*lTX»XITAB»TETAB»jXTAB.CTX,lXP0S> 
call xntErp t x+uxrho »sm i .smIp » i sx , x i tab ,smtab , i xtab ,csa , i xpos) 
jU ■ 15H5.0/SMI 

Z2 • ZH/ZS/Ri 

2C7 RhoEp*(Z2-Z1 >/UXKH 0 *Z3 

2 l 0 RHOUE » rhoe*ue 

KHOUCP“RHOE*UEP*OE*RhOEP 

ZHO“ZhOO*(TE/TG)**Z^VlS 

HO * HE ♦ 0L2U2 ....... . ... 

HOP ■ HEP ♦ UL • UE P 
PRE103 ■ PRAf 4 0T**I 1*0/3. 0) 

HAW»HE*PREl03*UE202 

CALL SEVAL12,TaW,E k ASE3.HAW) 

IFtlTATAB) 11.12.13 

11 T W® T A w 

H >» ■ H A A 

HWP-HEP+PRE 103*UE«uEP 
CO TO 1 H 

12 Tn®T W f AB ( l ) 

MA*T*TAB ( 2) 

HsiP*0 • . 

CO TO In 

13 CALL XNTERP ( X, U. TWP, ITWX, X 1 T Ab . TwTaB, 1XTAB. CT«X, 1XP0S ) 
CALX SEVALI V.Tm.CPw . HA V 

HWP»FJG*CPW*TWP 

f H IF ( T w .LE. TAW) "<j'0 TO 1 70 

riKlTE16.1S5» TA.TAft ...... .... 

IBS FORMAT < HbrtO**BARpRo pA "l LURE . « * WALL TEMPERATURE IT**. FB.2. 

1 62H ) CALCOLATtO <;REATER THAN ADIABATIC WALL TEMPERATURE ( TAW ■ 

2 i F 8 » 2 . 3 h T7 i 

WRl TE 1 6 . 1 06 I X, ZME» THETA, PHI 

nr<T FORMAT 1 23h AX 1 aL DISTANCE X = , IPEIH.7, SX , IIhMaCH NO. ■ , 

1 E1H, 7, SX, BHTHETa 1 EIH.7 . 5X, 6 H P H I 1 - . E 1 H . 7 I 

WRITE 16,2501 

250 FORMAT ( 6HH .CHEqH CONTOUR ANO MACH NUMBER DISTRIBUTION TAuLES F 
1 OR ERRORS. / UCh *MORE INPUT POINTS MAT BE REQUIRED TO ADEWUATEL 
2 Y DESCRIBE DERIVATjVE VALUES ALONG THE CONTOyR AT THIS POINT. / 
396H * A SMALLER RU n Ge-KUTTA STEP SIZE MAT BE REQUIRED TO ADEhUATEL 

H Y approximate integration values. // > 

CALL «Ul TS 

1 70 A « H* 

B»HO-HW 
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C»-UE202 
TFINJ-TE 
3 CALL ZETAlT 

crey=khque/zmu 

RThE»CREY*THE TA 

RPh i ■cr£y*ph“1 

CF = u»02bQ/(KTHE**0.2Sa) 

CF A G ■ 0.O2&U/(RPHi**O.2SO) 

CnPARl ■ l.C - PKA N OT ♦ A LOG ( 6*0 / ( &.J*PRANOT ♦ l.UJ) 
Crt*(PHl/TH£TA> **ZRsTaN*(CFAG/ 2.» / ( J •»5.*Sh|RT(CFAG/2» )»CHPARl) 

IF tlTATAB .LT. 0 > CH ■ Q.O 
£ R A S E l a R HOUEP/R M (J UX 
ERAS£2* < 1 .♦OElSOT ) /UE*U tP 

CALL XNTErp ( A, Yr, YRP, IYa, XITAB, YITAB, IX7AD, CYA, 1XFUS ) 
DARCsSQRT ( 1 .■.YRPaYftP) 

COF0RC*RH0UE/G*0E/dARC*CF/2. 

IF (EPSZ .LE. C.C ) GO 10 HQ 
EFa SeF i' ERA S E 1 ♦ E / Y K * Y R P 
Hu THETAP * CF / 2 • 0*0 ArC - T HE T A * ( EK A S E 2 ♦ ER ASE 1 ) 

EKASE2-H0-H* 

PlilP « CM*0 ARC/ERAsE 2* (HAW-HH ) - PHl*(£RASEl - ( H0P-HRP ) / E R A S E 2 ) 

i f ( i ho • he • n rlTukh 

IF (ilftTAo .LT. £) GO TO 66 
"Xi¥'"i - XHW£7rj»Cri/G* ( HA*”- "HA ) ' 

Hu«WA/(TAA-T») 

66 yOAO ■ GO A 

OFORCO-PFORCE 
OFLATu«OFlaT 

IF (EPSZ .LE. 0.0) GO TO HH 

ERA SET ■ P IT*Tr 

GOA ■ ERAsE 1 • G*< 

uFORct = ER ASE Y* C OF Jr C” 

OF L A T »0 . 

GO To HS 
H H GOA a 0 *• 

OFo R CE «C~d7 0 K C / 2 ~» 
oflat=oforce*ykp 
h5 YOARC » Y 2 ARC 
Y2ARC*0ARC 

IF <0A .LE. 0.0) RETURN 

Call XNTERPU - OX/2.0, EK ASt l , ERASE2, 1 Y X , XITab, YITAB, 

A 1 A T A B i CYA, IXPOS) 

YIARCaSGRI(l.+ERAS E 2*EKASE2) 

DXLARCa( YuARC + H.*Y j ARC-.Y2ARC ) /6.*0X 

XLARC=XLARC+OXLAKC 

SUMGOAaSUMG D A+OXLARC*(yOA+GOAO) 

IF (ICOOL .EG. 0) GO TO 2 
SO MGGA « COEFCL*SU m Q'oa ' 

CALL XNTERP l T l l . Z m Yul ,ZP , IZx ,TZT AB » ZmYTab , I 1 ZT AB ,CZx , 1 TP0S> 
1TP0S * IZX 

OlATUB ■ 2.C»SGKT( At/PlU 

REYL a MASSL*OIATU B /iAL*TOBEN*ZMYUL) 

J_ FORCE * __F ORCF ♦ DXlA«C*JOFORCE + OFORCO) 

F l A T a F L A T ♦OX L A K C * ( o F l A T ♦ D F L A T 0 ) 

RETURN 

S RALN a CREY*XLARC 
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ROLS=CREY*OELSTK 
IF (ZETA .GE. X .0 ) GO TO 62 
1 ■ l 
Z 1 *Z I 4 
Z2« Zlb 
Z 3 = Z I 6 
24»ZI 7 
Z5-ZI IP 
GO TO 63 
62 1=6 
Zl*Zl 1 

ZTSIT2 

Z3=Z I 3 
Z 4 ■ Z I 2 P 
Z5-ZI 3P 

63 A R I TE ( 6 , 5 1 ) 

5) FORMAT ( 1 HO/ 1 X , l 8 h c 0nT0UK PROPERTIES. 5A,15hFlOW PkOPeRTIES.OX, 

A THh^OUNDARy LAYER. VX, 13HHEAT TRAnSFEw»7X, 

X IBriK'llERNAL INTEGRALS, 7X,12HCQEFFlCIfcNTb/) 

WRITE I 6. 52) X.ZMt.oEulA.HG.ZETA.CF 

52 Format (ix,7ha = ,fi i.6,3x»7hzme ■,fj2.6,3a,7hoelta =, 

A 1 PE 1 3 > 6 » 3X f 7hHG ■ , OPF I 2 . 6 , JX , 6HZE 1 A » , I PE 1 4 . 6 , 3X . 

X 6HC.E =, 1 p E 1 3 » 6 ) _ 

WRITE (6 .53) XeARC.tE.BDELTA,G6*ZINTPRU ) .21 ,CH 

53 FORMAT (IX. 7 HXLARC ■ , f I 1 . 6 . 3X , 7H T£ • , F l 2 • 6 , 3 X , 7HB0EL T A« , 

A 1 P E 13 • o » 3 A , 7 Hi» * •, 1PE12.6.3X.A6, IkE14.6,3x.6HCH «. 

X 1 PE 1 3 • 6 ) 

WRITE16.54) YR,TW, u ELSTK,SUMgOA,ZINTPk<1«-D,22,KTH£ 

5 4 FORMAT (jX./HYK »,F l 1.6.3X.7HTW = , F 1 2 • 6 . 3X , /HDELS T R= , 

A"' ' J PE I 3.6.3 A , 7 h SUM 90 A" .IPE12.6.3X.A6, 1 PE l 4 . 6 , 3X , 6HRTHE =, 

X IPE13.6) 

WRITE (6, 5b) YRP,TA*,TH£TA,FORCE,ZInTPR< l ♦ 2 ) ,t3,RXLN 
bb FORMAT <Ia.7hYKP ■ . F I I . 7 . 3 X . 7 H T A A ■ » F l 2 • 6 , 3 X , 7 hT HE T A =. 

A lPE13.6.3X,7hFQKCE * , OP? 1 2 . 6 , 3A . A 6 , I r E 1 4 . 6 , 3 X . 6HR XLN =, 

X IPE13.6) 

YKUELb * YR - 0 E L b t R 

WRITE (6,561 YKOtLs .ZhEP , PHI , FUAT , Z INTPRI 1*3) ,Z4 ,KPHI 

56 FORMAT I IX,7hYhUtLs«,F H.e,3X,7HZMEP ■ , F 1 2 . a , 3 * , 7hPH 1 », 

l l PE 1 3 • 6 , 3 X • 7 HKL aT . . CPF l 2 • 6 . 3X t A6 . 1 PE 1 4 . 6 • 3X • 6hKPH I ..1PE13.6) 
WRITE (6,57) UE,DE L SoT , RB A R . Z l N T P R ( 1*4) .Zb.RwLS 

57 FORMAT (22X.7HUE » , F 1 2 . 6 , 3 X , 7H0 E ESO T ■ . F l 3 . 6 , 3* , 7HRb aK »»FI2»6, 

A ~ 3x7 A 6 . I PE ) 4 V 6 , 3 X , 6HRDUS =,1PE13.6) 

WRITE (6,58) Pt.RhoE.PRANOT.GAME.SMOL 

58 FORMAT ( 2 2 X » 7 h pE = , 1 pE 1 2 • 6 , 3A . 7HKH0E * . CpF l 3 . 6 , 3x . 7 HPR A NO T» . 

A Fi2.10.3X.6Hfc.AME «,FI4,8.3X,6 hSM0L « ,F 1 3 .6 ) 

CUSAL * I.O/OARt 

IT (EpSZ .lE. 0*0 ou TO SOO 

OLLF 1 ■ 2.0*PlE*YK.RhOUE*ThETA*0E*C0bAU/G 

OtLF 2 ■ 1.0 - oELS u T*P£/ ( KHOOEAO t/G > 

OtLFA ■ 0ELF1*0ELF 2 

THRUST ■ PIE*YR*«2* (RHOUE*UE/G + PE » 

0 EFT hR * 10£s«0*UELpA/THRUbT 

ThRUSA ■ PlE*YR*»2 * ( RHoUE*UE/G ♦ PE - 2)16.2240) 

ZHAbSR ■ P)E*(YK ♦ OtLSTR.COSAL ) •»2*HH0UE 
XMAbSR ■ PIt*YK»*2*RH0UE 
GU TO 510 
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C THE TWO OlHEWSnWXE - c*5E'~ ASSUMES A WIDTH OF ONt FOOT 

C 

T>GC OELF1 * 2 • Q* RFiQUE* t RET A*UE*COSAL/<i 

DELF2 * 1.0 “ UELS 0 T*PE/1RH0UE*UE/U) 

DETLFa « “DELF l*OtLF 2 

THfiUSl » 2.0«TK* 1 HhOuE*UE/G ♦ PE) 

OEFT HR • 100.0*DEL>A/THK us t 

THRUSA » 2.0*TK* IRh 0uE *UE/G «■ PE - 21)6. 2ZH0) 

ZMASSK ■ 2 . 0* ( TR ♦ D E L S T R * C 0 S A L > * R H 0 u E 
XMASSK ■ 2.0»VK*RH 0 Ut 
5T5 VlSP ■ THROST/ZMASsR 
AISP * THKUSa/ZMASsR 
TbLlSP ■ -DEEFa/TLqCkT 

WRITE 16,1) AISP. X M ASSR, ThRUSA. OELFA,T bLiSP,vlSP,ZMASSR.THR u ST» 

1 ' DEPTHS, CuSAL 

1 FORMAT I//2bX ,blHT H RUST DEFICIENCY AND SPECIrIC IMPULSE DECREMENT 
t UUE~TO THE BOUNDARY E A Y E R £ FT'ecTTTs A 7 6 HA ISP F 1 0 . H.S X , BH A M AS SR * , 

2 F 13. H.5X.BH THrUSA *»F 1R»H»&X»7H0ElF A ■ »F l 2 . M , S A . bHTBL 1 SP »»F11.6/ 

3 sX, 6 hVISP ■ ,F fO • S , 5 A , 8HZM ASSR «,F1 J. R , 5 A , aH I HRUS T *,F1H.*»,SX, 

1 8H0EFTHR ■.F11.8isX»7HC0SAL »,F12.8//) 

if ticoot .Ey. o> Return 

TwL a TL 1 

Call anterp (TlI.CpE.cpp , i c x , tztab , cpl t ab , i tztab , c cx ,i t pos ) 

CAUL ANTERP <TL1 ,K a M0L,KP,1KX,TZTAB.KAMTAu.I IZTau.CRX.ITPOS) 

CALL ANTXrP r TLTTZ hY UL, Z P ,1 Zx iTZTAB » ZmYTAB , I |ZTAb,CZx, 1 TPOSl 
CALL XNTERP <X,ENrt A ,tNHAP,l£X,XITAB.tNHTAB,lATAU,CEX,l *POS> 

pkandl ■ cpl.zmyul/Ramdu 

7 C TAL6 » TWL 

HL * 0 *02 s a * ttaru L / U I A 1 0 3^ R E TL • • 0 • B 0 * P K A nD l * • u • 0 • l TL 1 / TWL) ••0.550 

HL a M L * E N H A 

SA 1 * HL*U»C RAhDw/ (THICK*H$) ) 

S A 2 a RAMD»v/TMICK 

T«L » (Sa1*Tli ♦ S A 2*TA« ) / ( SA 1 + 5A2) 

IF ( AbS I T WLO - TWL) • (> T « 0.010) GO TO 70 
TEmPrl * TWL/Tcl 

TrtbCA ■ <HG*TA« ♦ rAMDW/THI CK*T»L ) / IHG ♦ KaMuW/THICK) 

U.1 1 a HSMTAw - T A(, c A > 

IF IEPSZ .Lt. 0.0) GO TO 600 

SST a PI£.YR*DWI *U E L Ab A #D ARC aCOEFCL 

«D To 6 10 

6oD SST a DwT*DtLABA*D A RC»COEFCL 

610 TLCA a ( TL l ♦ T L 2 ) / 2 • Q «■ SS T / I CPL*M A SSL ) 

IF ( I DUMP ,GT. 0) TLCA ■ (TLO + TLD/2.0 + SST/ (CPL*M A SSL ) 

OMASSi. ■ MASSL/ ( A'L*TUbEN) 

SUMO* 1 * SUMO* I ♦ 1*2.0 

WRITE 16,711 DmASSl»HL.DWI , RE Y L . SUMO* l . TE mPRu , TLC A , TwGCA , T wL . 

1 ' ' UlATUg, THICK, SUMDGA 

71 FORMAT </SQX,31HRE 6 EnERATIVE CALCULATION OU Tp U T/ / S X , OHO* A SSL a, 
l F12*M,SX»HHHL *»F)0»6,7x»5Hii|Wl » . F 1 2 . R , 5 a , 6mRE YL «.F19«H.SX, 
z BhSumOWI *tF lS.6/sX,bHTEMPRL * »F 1 0 « M , 7 A , 6H TlC A » , F 1 o* <t i SX , 

3 7HT wuC A a,F lCt.H,&A»SHTWL a , F 1 0 . 1 , 1 SX , 8HD 1 A Tub »,FIS,10/ 

1 SX , 7hTH I CK =.FIG» 6»B1X i8HSUMPGA ■ , F 1 S . 6// / / / / / / ) 

RETURN 

END 
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barset 


subroutine barset 

COMMON /CSEVaL/ NOC T AB i I S . H 0 J . F J 6 , C J 6 , GM I 02 , GOGm 1 , POM A X , C P 0 . HO • 
A SO,TCTAa(2CI.CPTAB(2Q).aCP<20).CCP(Z0).0CP(2a>. 

K o T A B ( 2C ) , H T A B l 2 0 ) ,BARBl(20) ,BARn2(20) ,BaRB3(20) 

COMMON /INPUT/ I 0 A M A X , I C T A B • I P R I N T , I T A T A B i I X T A B . M ZE T A . D *h A X , 

A EPSZ.Fj.G.GAMO.PQ.PHll, PIE. PKANOT.KB aR. SCALE, TO* 

K THETA l ,TOLCFA.TOlZET,TOLZME,ZMUO.ZMV15,ZNSTAN 

common /Tables/ PtTAB! 100) .SMTAB! ICO) ,TETAB< 100) .THTAB! 100) . 

1 OETABl loo) .XITABI 100) ,YITAB( 100) .ZMTAB1 100) 

pie«3. 1415926s 
F J6«F J »G 
ROJ.RBAK/F J 
TMAX.TC 
I 1" l 

IF(ITATAB) 15,12.11 
tl U«IXTAB 
12 00 14 1.1,11 

IF UwTABtl) .LE. TMAA) GO TO 14 
TMAX . T.1TAB1I) 

14 CONTINUE 

15 IF I ICTab .Ed. 0) 00 TO 20 

IF I T M A X .LT. TCTAB t ICTAB ) ) SO TO 16 
<*RITE(6,52) TMAX , TCTAB t ICTAB ) 

52 FOHMAT ( // 44H **BAR5ET ERROR** TEMPERATURE InPUT VALUE (, 

1 1PE14.7, 29H) EXCEEDS TABLE UPPER LIMIT (, E14.7, 1H) ) 

60 TO 57 
16 NOCTaB*ICTAB 

CALL BMF ITS ( TC T AB , CPT AB . I C T AB . BCP , CCP . OCP ) 

TIEl. TCTAB! 1 ) 

TIE2«TIE1*TIE1 
TIE3-T !E2*T1E1 

HTAB(I) . CPTAB(1)*TIE1 - BCP1 1 )*TlE2/2,0 

BARaKl) « CPT AB ! 1 • - BCP(1)*TIEI ♦ CCP(1>*T1E2 - 0CP(1)*TIE3 
BAKB2 ( 1 ) ■ B C p ( 1 ) - 2 . C*CC P ( 1 > • T I E 1 ♦ 3 . 0*DC P I 1 ) • T I E 2 
BARB!! I l*(CCPt 1 ) - 3 • .DC P ( 1 ) • T I E 1 )/2, 

stab! 1 »*-BARBl < 1 )»AL0G!TIE1)-BARB2< 1 ) * T I E 1 -8 AR83 ( 1 >*TIE2 
X -OCP! ) )/3.*TlE3 

6)>0. 

00 14 I«2. ICTAB 
THE J»TlEJ 
THE 2 .T I E 2 
TmE3*T I E 3 
TlEl.TCTAB! I I 
TIE2.TIE1.TIEj 
T 1E3.TIE2.TIE1 
OELT«TlEl-TMEl 

HTAB(I) . HTABtl-l) ♦ CPTAB ( I -l ) .CELT ♦ BCP ( I . 1 » .DEL T . . 2 / 2 . 0 ♦ 
A CCP< 1-1 )*OElT*. 3/3.0 ♦ DCP( I"1 )*OELT». 4/4.0 

IF (1 ,GE. ICTAB) GO TO 19 

BARBl(l) ■ CPTABIl) - BCP ( 1 ) *T I E 1 ♦ CCP(l)»TIE2 . 0CP(I|.T1E3 
8ARB2 ( 1 ) «BCP( I I -CCP ( I ) / .5»T IE 1*3.»0CP ( 1 ) »T IE2 


BARS 1 

/CSEVAL/ 

/CSEVAL/ 

/CSEVAL/ 

/INPUT/ 

/Input/ 

/ i nput / 

/TABLES/ 

/TABLES/ 

BARS 24 
BARS 25 
BARS 26 
BARS 27 
BARS 2B 
BARS 29 
BARS 30 
BARS 31 


BARS 34 


BARS 39 
BARS 40 
BARS 41 
BARS 42 

BARS 51 
BARS 52 
BARS 53 


BARS 57 
BARS 58 
BARS 59 
BARS 60 
BARS 61 
BARS 62 
BARS 63 
BARS 64 
BARS 65 
BARS 66 
BARS 67 
BARS 68 


BARS 73 
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BaRB 3 ( I 1-lCCPl lJ-3**DCPil I •TIEJ )/2. 

G1-gY*BARB 1( 1-1 I*AL0G(T1E1/TME1 >*BARB2t l-ll*DELT 
X *BARB3(I*11*(TIE2-TNE21*DCP( 1-11/3«*(TIE3’TME31 

G 2 »BARBl (I )*AL06(TIE||*BARB2(l l*TIEl*BARB3t I >*TlE2*OCP(ll/3.*TlE3 

GT AB (I )a«l*G2 
19 CONTINUE 
IS-ICTAB-l 

CALL SEVALI I .TO.CPO.HOl 
SAMOaCPO/ (CPQ-ROJ) 


BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 


20 IF (GAMO *GT* I • 0 » GO TO 56 
WRITE (6.591 GAHO 
FORMAT ( Bl* **BaRSET ERROR** 

1 greater than one C 1 » • gamq ■ • 

2 UNITS**'* CP* RBAR, FJ // > 

S7 CALL QUITS 

56 GM102 ■ (GAHO - 1.0)/2*0 

GOGH l *G AMO/ (GAHO- l • > 

POMAX»PO 

IF (1CTAB *GT, 0) GO To 30 

NOCTAB « 6 

NOCTM1-NOCTAB-1 

I S*NOCTM 1 

CPQ-GOGMl/F J«RBAR 

CJG«CP0*FJG 

H0*CJG*T0 

00 23 1 * I , 1 XT AB 


RATIO OF SPECIFIC HEATS MUST BE BARS 
E 1 9 *7 / 96H CHECK FOR I NCONS l STENTB ARS 

BARS 

BARS 

BARS 

BARS 


BARS 

BARS 

BARS 

BARS 

BARS 

BARS 


TE * TETABI 11 

IF (TE *LE. TMAXl GO TO 23 
TMAX • TE 
23 continue 

TCTAB(NOCTAB)*TMAX*100. 

TCTABC 1 »*1 *e-io 
Z l-NOCTMl 

OELT - (TCTAB (NOCTAB )-TCTAB (11 )/Z l 
ERASE 1*-CP0*AL0G(TCTAB! 1 > > 

DO 25 lal.NOCTAB 
GT AB ( I taERASEt 
CPTABI I I * CPO 
BCP( l ) ■ 0 * 

CCP( 1 )aO. 

DCP ( 1 )*0« 

BARB 1 < I l*CPO 
BARB2 ( 1 )aO* 

BARB3 ( 1 )*0* 

HTAB ( I)aCPO*(TCTAB( I )-TCTAB( 1 ) ) 
IF (I • GE « N0CTM1I 60 TO 25 
TCTAB(I*1» ■ TCTAB(l) ♦ DELT 
25 continue 

30 IF (ITWTAB ,NE. 01 GO TO 38 

call sevalii. twtabiii* erasei* 


BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 

BARS 


BARS 


TWTAB12) ) 


38 CALL SEvALlO*TO*PO*SO) 
RETURN 
END 


BARS 

BARS 


7 ** 
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76 

77 

78 

79 

80 
81 
82 


86 

87 

88 
89 

92 

93 


99 

100 

101 

102 

103 

109 


108 

109 

110 
111 
112 
113 

119 
115 

117 

118 
1 19 

120 
121 
122 
123 


126 


132 

133 
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bmfits 

SUBROUTINE bmfits I X . Y • N • BL • CL • DL I 

DIMENSION M20> .B(20> ,C(20> ,FC20> ,GI20! ,X(2Q> ,Y<20> , 

A BH20I .CL I 20 .DL(20> .FL120I , YPP 1 20) 

C 

1 » 1 

It FLU) * XlM) - XII) 

1*1 + 1 

IF tl .LT. N) GO TO 11 
1*2 

IS BUI • -FL( 1-1 >/FL( I ) 

All) » -2 • 0* ( FL ( 1 ) ♦ FLC l-l > )/FL( I) 

C(l! ■ 6.0/FLl I 1*1 I Y< l + l I - Ytlll/FLIll - ( Y ( I ) - Y I I - 1 ) » /FL I I - 1 1 ) 
1*1 + 1 

IF II .LT. N) GO TO IS 
G ( 2 ) « 1.0 
F (21 ■ 0.0 
I • 3 

02 GUI ■ All'll ♦ B( 1-ll/GI I-n 

F(I) ■ -IBII-l » *F i I-l)/Gll-t) + CII-DI 
1*1 + 1 

IF II .Lt. N) GO TO 32 
YPPINJ • FINI/IGIN) - 1.01 
YPP(N-I) ■ YPPIN1 
I ■ N - 2 

M7 YPP I I I « (YPPIIMI + F I I + l II /G( 1 + 1 I 
1 * 1-1 

IF II .GT. 0) GO TO H7 
I * 1 

SI Be I I I » I Y I I ♦ 1 • -Y I 1 1 I /Ft< I I - IFL « I I • I YPP I I + 1 I ♦ 2 • 0* YPP C I ) I » /* • 0 
CLI I 1 » YPP (I 1/2.0 

OLll) ■ IYPPU + 11 - YPPI I 1 I/I6.0*FLI1) 1 
1*1 + 1 

IF I 1 .LT • N) GO TO St 

RETURN 

end 


78 



>J (MT X Oi N 


cflval 


CFEV 


C 


C 


function cfeval <cf ht i 


i 


DIMENSION X<8)»A<7),8(7),C17I«0<7),IX(8| CFEV 2 

EQUIVALENCE (X. IX) .12.12) CFE V 3 

DATA A/-2. 079177 3E-2.-4.971S425E -3, 1 . 26 l 4 39 2E- 3 . - 1 • 008 8 6 I 7 E- 3 . 

1 1.752I422E-4.-2.883630E-4.5.998S794E-6/, B/ 0 • 209 1 586 2 , 

7.3894S60E-2.-I.6799227E-2.2.95l99nE-2.-1.862082lE-3. 
I.33l87H7E-2»2»lC3b7o7E-3/» C/*0»92l4204 3. -0*535 92107, 

-9. 607S3CE -2.-0. Ml 101 | 15. -0.1 390*1416. “0.29826897. 

-C» 15583627 / 1 D/-M.MM57710.-M.8 l 19952 .-5.5230787. -4. 8 09 2254. 
-5.6024598,-5.0345585.-5.63752327. J/17. x/2. 5099998. 

1 7. 287782, 127. 74039, 897. 8M729. 6310. 6880. 44355. 857, 

8 327747 . 91 , 1982759 . 2 /, ZERo/O. 0 / 


2-CFRT 

CFEV 

IS 

IF (2 .LE. 0.0) 50 TO 3 

1 IF < I 2-1 XI JU 2,7.9 

CFEV 

17 

2 J« J- 1 

CFEV 

18 

IF t J ) 3, 5, 1 

CFEV 

19 

3 J» 1 

CFEV 

20 

6RlTt(6,4) Z , ZERO , X ( 8 ) 

CFEV 

21 

4 FORMAT 1 / 1 OX , 1 4HCFEV AL F A I LURE » 5 X , 3HZ 
lR0M,5X,E18.8.2X,2hT0,El8.8) 

CALL QUITS 

», 1PE1S.8.5X. lbHLlMlTS A*E F 

CFEV 

23 

5 IF 12 .Lb. 0.0) 50 TO 3 

J » 1 

Y«.009896/Z»».562 

CFEV 

26 

GO TO 8 

CFEV 

27 

7 ZL»AL0G(2) 

CFEV 

28 

YL-0(j)+2L*(C<J)+ZL*(B(J)+ZL*A(J))) 

CFEV 

29 

Y«EXP(YL) 

CFEV 

30 

8 CFEVAL-Y 

CFEV 

31 

RETURN 

CFEV 

32 

9 IF 1 IZ .LE • IXIJ+1 ) ) GO TO 7 

J ■ j ♦ 1 

IF(J-B) 9,3,3 

CFEV 

3S 

END 

CFEV 

36 
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DIRECT 


C 


SUBROUTINE DIRECT 

10 CALL REaDIN 
CALL BARSET 
CALL BARCON 



RE l 


RE 2 
RE 3 
RE ** 

RE 7 



FI IF 

FUNCTION FIIF (SI 

C 

COMMON /COFIIF/ I F I NT , AF l N T , B F I NT , CF I NT , MM l NT , T f I NT 

STOM"! . 

IF (MMINT •«.£• 01 GO TO 12 

00 *» I«l, MMINT 

‘I STQM ■ SI.OM* S 

12 FDEn ■ AFINT ♦ S«IBFINT ♦ S*CFInT) 

1 E_ .1 LF. I ALT * G-E-t ..ZJ <LU_T0 2 

FNUM ■ STOM* I t.o > SI 

60 TO a 

2 FNUMmSTOM 

3 CAUL SEVAII2.TtO.FOEN) 

FI IF-FNUM/T*TF1NT 
Return 

end 


/COFIIF/ 

FIIF ** 

FIIF 7 
FIIF 8 


FIIF 13 
FIIF IH 
FIIF IS 
FIIF 18 
FIIF 17 
FIIF 18 
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GETPT 


subroutine getpt (zme.pi.tii 

COMMON /CSEVAL/ NOC T AB . 1 S , RO J , F JG , C JG , GM l 02 , GOGM 1 , POM A X , CPO , HO » 
A SO i TCTAB ( 20 ) > CPT AB ( 20 ) . BCP < 20 I . CCP ( 20 I , OCP ( 20 ) * 

K GTABI20). HTA8< 20 ), BARBU 20 1, BARB 2(20), BARB 3( 20 > 

COMMON /INPUT/ IOXMAX.IcTaB.IPRINT.ITWTaB.'IXTAB.MZETA.DXMaX. 

A EPSZ.Fj.G.GAMO.PO.PHlI . P I E . PR ANDT , RB AR . Sc ALE > TO • 

K THETa I .TOLCFA.TOLZET.TOLZME.ZMUO.ZMVIS.ZNSTAN 

ZME2«ZME»ZME 
PR001-2./RBAR/ZME2/G 
0ENM2* 1 .♦GMI02*ZME2 
TE«TQ/DENM2 

IF (ICTA8 .GT. 0) GO TO 20 
PE*PQ/DENM2**G0GMi 
IS PI * PE 
T l«TE 
RETURN 

20 ITER m 0 

TOL ■ TOLZME/ZME 
21 TEO-TEG 
TCO«TC 
TEG»T£ 

CALL SEVALI 1,TE,CPE,H£> 

GaME«CPE/ ICPE-ROJ) 

TCMHO-HE)/GAME*PRODI 

IF (ABSMTC - TEl/TE) .LE. T0L1 GO TO 30 
IF (ITER .GT. 0) GO TO 2H 
TE • I 2.0*TE ♦ TC1/3.0 
GO TO Z8 

2 H IF (ITER .LE. SO) GO To 27 

RRITE ( 6 . 26 ) ZHE.TC.TcO.TE.TEO 

26 format ( 31ho**getpt failure... mach no. - . ipeih.7 / ihx. 

1 1 7HT (CALCULATED) * , 2E16.7 / IHX, 17hT (GUESSED) 

2 2EI6.7 // ) 

GO TO 30 

27 ZK»( TC-TCOI/ITE-TEO) 

TE»(TC-ZK*TE) / ( l .-ZK ) 

IF (ABSIITE - TEG ) /TE ) .LT. TOL) GO TO 2? 

IF (ITER .LT. 10) GO To 2S 

IF ( ABS | ( TE - TEO)/TEi .LT. TOL) GO TO 29 

28 ITER«ITER*l 
GO TO 2 1 

29 TE«(TE*TEG>/2. 

30 CALL SEVALI-I ,TE .PE. 50) 

GO TO IS 

end 


/CSEVAL/ 

/CSEVAL/ 

/CSEVAL/ 

/INPUT/ 

/INPUT/ 

/INPUT/ 

GETP 9 
6 ETP 10 
GETP II 
GETP 12 

GETP IS 

GETP 18 
6ETP 19 


GETP 23 
GETP 28 
GETP 2S 
GETP 26 
GETP 27 
GETP 28 


GETP 32 


GETP 35 
.GETP 36 
GETP 37 
GETP 38 
GETP 39 
GETP HO 


GETP HH 
GETP HS 
GETP H6 

GETP H? 
GETP SO 
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INTZET 


SUBROUTINE INTZET (X1.X2.ZINT) 


- 

INTZET - QUADRATIC FOUR POINT I NTEGR AT I ON SCHEME 

- ZETA(BARTZ)INTZ 

2 


DIMENSION XC(2»> ,YC(2l ) ,Ym(4) 

INTZ 

3 


DX2 I -X2-X 1 

INTZ 

7 


SUMINT-O. 

IF ( DX2 1 .Eq. 0.0) 60 TO IS 

INTZ 

8 


DXC-OX21/20. 

INTZ 

1 1 


I MAX--9999 

INTZ 

12 


FmAX»-1 .E30 

INTZ 

13 


DO 10 1-1,21 

INTZ 

IM 


XC( 1 l-Xl+FLOATC 1-1 ) *DXC 

INTZ 

15 


YC ( I ) -F i I F ( XC ( I ) ) 

IF ( YC (I) .LE. FMAX ) GO TO 10 
IMAX > I 

INTZ 

16 


XmAX-XC t 1 ) 

INTZ 

19 


FMAX-YC ( I ) 

INTZ 

20 

10 

CONTINUE 

IF (0X21 .GT. 0.10) GO TO 17 

INTZ 

22 


SUMINT«10.»YC(l)»16.*YC(2)-2.»YC13) 

INTZ 

25 


DO 14 1-2,19 

INTZ 

26 

IM 

PaRINT*13.»(YC<1)»YC(I*1))-YC(1-1)-YC(1+2) 
SUMINT ■ SUMINT ♦ PARINT 

INTZ 

27 


SUMINT-SUMlNT+10.»YC(2lj>16.»YC(20)-2,*YC (19) 

INTZ 

30 


SUMINT-SUMINT/24.-DXC 

INTZ 

31 

\s 

ZINT-SUmINT 

INTZ 

32 

17 

return 

FB«K ■ FMAX-0.20 

INTZ 

33 


SUAlNT*0. 

INTZ 

35 


subint-o. 

IF UMAX .LE. 2) GO TO 21 

INTZ 

36 


DO 19 1. 2. IMAX 

IF ( YC ( I ) .LE. FBRK ) GO TO 19 
IBRK * I - 1 

INTZ 

39 


60 TO 20 

INTZ 

M2 

1 9 

CONTINUE 

INTZ 

M3 

20 

ibkk-imax-i 

IF ( IBRk .GT. 1 ) GO TO 22 

INTZ 

MM 

21 

I BRK- 1 

INTZ 

M6 


I BRKM 1 >0 

INTZ 

M7 


60 TO 2S 

INTZ 

M8 

22 

IBRKMI-IBRK-I 

INTZ 

M9 


SUAlNT“I0.»YC(l)+16.*yC(2)-2.*YC(3) 
IF (IBRK ,L£, 2) GO TO 204 

INTZ 

SO 


DO 23 I.2.IBRKMI 

INTZ 

53 

23 

20M 

25 

PAR I NT" 13.* (YCl I )+YC( I ♦ 1 1 )-YCl 1-1 )-YC( I*2> 
SUA1NT ■ SUAINT ♦ PARINT 
SUAINT ■ SUAINT/24.C*DXC 
DXM » DXC/3.0 

IF (IBRKM1 .GT. 0) GO TO 206 
K » 2 

INTZ 

5M 


JS ■ 2 

INTZ 

63 


60 TO 207 

INTZ 

6M 
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206 K - 3 
JS • 1 

07 00 26 I » 2»H 

X M*X C ( IBKK)+FL0AT(1-K)*0XM 

26 YHIH ■ FIIFCXHI 

IF (IBRKMl • G T • 01 GO TO 209 

SUB 1 NT ■ I0.0*YM12) ♦ U.0*YM<3) - 2.0*YM(H) 
09 DO 27 I * I BHK t 1 9 
DO 28 J * J S * 3 
YMt l )«YM12> 

Ym(2)*YmI3) 

YM(3)»YM(M) 

XM»XM + D X M 
YM(H)«FI IF <XM) 

HARINT«13.*(YM<2)+YM3>)-YM(1)-YM(<4) 

28 SUBINT . SUBINT ♦ PARlNT 
JS » 1 

27 XM * XC ( I + 1 ) ♦ D XM 

DO 2V J»1 .2 

YM( 1 i ■ Y M ( 2 I 
Y M ( 2 ) ■ YM ( 3 ) 

YM(3)«YM(H) 

XM*XM+OXM 
YHCHImFI IF tXM) 

PARlNT»13,«(YMl2I+YM(3))-Yhm-YMm 

29 SUBINT « SUBINT ♦ PARlNT 
SU8lNT»SUBlNT+10**YMlH)**-I6.»YM(3)-2»*YMt2) 
SU61NT«SUB INT/2*U#DXM 
SUMInT«5UAINT^SUB1NT 

60 TO IS 
END 


hajntb - - ...... 

C I C R P 6 REFERENCE PROGRAM TBL 

C DECK SEOUENCEO BY SUBROUTINE 

C 

COMMON ./INPUT/ IOXMaX*ICTaB*IPR1NT,IT«TaB»IXTaB.M2ETa.DXMaX, /input/ 


a epsz.fj.g.gamo.po.phii .pieiPRanot.rbar, scale. to» /input/ 

K- . THETAl*TOLCFA,TOLZET,TOLZME,ZMUO.ZMVIS t ZNSTAN /INPUT/ 

c 

1.0XHAX a o 

CAUL DIRECT TBL I 

END TBL a 


85 



w -e u* K) 


quits 


QUIT 


1 


c 

c 


c 


c 


c 


c 


c 

c 


c 

c 


c 


subroutine quits 

COMMON /COFI1F/ IF 1NT , AFlNT .BFINT .CFINT .MMINT ,Tp INT 

COMMON /COOL/ I C OOL ■ I 0 UMP , I T Z T A B . AL . C OE F CL , CPL . OE L XB A . D I A TUB , 

I FLQ*RT.MASSL.PRANDL,RAMDL,RAMO#,R£YL.SUMQGA.SUMQl*I. 

THICK , TLO.TLI .TL2.TLCA .T0LITE.TUBEN.T6GCA .ZMYUL . 
CPLTAB ( 20 ) ,RAMTAB(20I ,TZTAB(20) ,ZmYTAB(20> , 

ALTABt 100) *THITAB( 100) .TLCTABI 100) tTLTABl 100) • 
T66IAB I 100) 

real massl 

COMMON /CSEVaL/ NOC T AB . I S , RO J , F JG , C JG , GM 1 0 2 , GO Gm 1 , POM A X , CP 0 , H 0 • 

A S0.TCTAB(20) .CPTAB120) • BCP ( 20 ) , CC P ( 20 ) , DCP ( 20 ) • 

K GTABI20) ,HTAB(20) .BARB l 120) >BARB2(20) ,B ARB 3 (20) 

COMMON /INPUT/ IOXMAX , ICTAB . IPRI NT. I TWTAB . IXTA8 .MZETA .OXMAX, 

A EPSZ.FJ.G.GAM0.P0.PH1 I .PlE.PRANDT.RBAR.SCALE.TO. 

K THETA I .TOLCFA.TOLZET.TOLZME.ZMUO.ZMVIS.ZNSTAN 

common /inter/ CFAGT.CFAGP.CHPARl .OX.DXRHO.HE.H*. ibeg.mzetam, 

A OOMZeT,PhIP.PRE103.RHOE.RhOUE,RMZETA,THETAP. 

K X l BASE, XI END, ZET ATM .zmzeta ,ZMz£Tm,ZMZETP 

common /LOOKUP/ ICX.ImX.IPX.iRX.ISX.ITPOS.ITAX.iTX.IUX.IXPOS.iTx, 

1 IZX.CCXI6) , C M X ( 6 ) « CP X ( 6 ) , C R X ( 6 ) < C S X ( 6 ) . C T A X ( 6 ) . 

2 CTXI6) .CUXI6) . C Y X ( 6 ) » C ZX ( 6 ) 

COMMON /NHANCE/ IEX.CEX16) .ENHTaBI 100) 

COMMON /OUTPUT/ BOELTa.cF.CH.OELTa.DELSOT.OELSTk, flat, FORCE. hg. 

A PE. PH I ,QI»,SUMQOA,TE.Th£TA.TA.UE,X,XlARC.YR,ZI . 

K Z2.Z3,ZR,ZS,ZETA,ZME 

COMMON /SAVED/ A.B.C.ZI1.Z11P.Z12.ZI2P.Z13.ZI3P.Z1R.Z1S.ZI6.Z17 

COMMON /TABLES/ PE T A B (I 00 ) • SM T A B ( I CO ) . T E T A B ( I 00 ) . T A T A B ( 1 00 ) . 

1 UETABI 100) .XITABI 100) .YITABI 100) .ZMTABl 100) 


WRITE (6. 1 ) 

l FORMAT ( 3MH l QU I Ts COMMON DIAGNOSTIC OUTPUT...) 

*R|TE(6,b) IFINT.AFINT.BFINT.CFINT.MMINT.TFINT 
b FORMAT ( //SOX .21HC0MM0N BLOCK /C OF 1 I F / / 25 X , U 0 , f P 3E 20 . 8 , I 1 0 , E 2<J • 8 
WRITE (6.2) IDXMaX, 1CTAB.IPR1NT.it A Tab. 1XTaB.MZeTa. 0XMAX.EPSZ.FJ, 
A G, GA MO ..PO.PHll, PIE, PRANOT.RBAR.se ALE, TO, T HE TAl, 

K TOLCFA .TOlZeT ,tolz«e .zmuo.zmvis.znstan 

2 FORMAT < //50X . 20HC0MM0N BLOCK / 1 NPU T / / / 3X • 6 I 6 • 5 < H X • 1 PE I 3 . 6 ) / SX . 

A 7«*tX.lPE13.6)/5x.7(MX.lPEl3.6)/) 

WRITE (B.lO) BOELTA.CF.CH, DELTA, OELSOT.OELSTR, FLAT, FORCE, HG, PE, 
a phi, qw.sumqda.te. Theta. tw,ue.x,xlakc.th,zi.Z2,Z3. 

K ZM.ZS.ZEtA.ZME 

jO FORMAT ( //SOX .21HC0MM0N BLOCK /OUT PUT / / / H ( &X , 7 < hX , 1 PE 1 3 . 6 ) / ) / ) 
WRITE (6,3) N0CTaB.1S,ROJ,FJG.CJG.GM102,gOGM1 .POMaX.CPO.HO.SO, 

A TC tab .CpTAB .BC p.CCp.OCp.GTAB ,HTAB ,B a RBI ,bARB2,BARB3 

3 FORMAT (//SOX.21HC0MMON BLOCK /C SE V AL / / / 3 X , 2 1 S , 1 P 9 E l 3 . 6 / 


/COFI IF/ 

/COOL/ 

/COOL/ 

/COOL/ 

/COOL/ 

/COOL/ 

/COOL/ 

/COOL/ 

/CSEVAL/ 

/CSEVAL/ 

/CSEVAL/ 

/INPUT/ 

/INPUT/ 

/INPUT/ 

/INTER/ 

/INTER/ 

/INTER/ 

/LOOKUP/ 

/LOOKUP/ 

/LOOKUP/ 

/NHANCE/ 

/OUTPUT/ 

/OUTPUT/ 

/OUTPUT/ 

/SAVEO/ 

/Tables/ 

/tables/ 

QUIT 7 
QUIT 8 
QUIT 9 

) 
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A < iPUXi 1PE12.6) ) ) 

WRITE (6,R) CFAGT.CFAGP.CHPARl , DX , D X RHO , HE , H W , I bEG , MZE T A M . OOM ZE T , 

A PHIP,PRE103,RH0E,RH0UE.RMZETa,THETAP,XIBASE. 

X XlENO .ZETaTm.ZMZETA.ZMZETM.ZMZETP 

H FORMAT J//SOX.20HC0MMON BLOCK / I NTER///SX ( 7 C 5X , jPE I 3. 6 ) /5X . 2 1 9 , 

A 6(5X«tPEl3.6)/5x,7(bX.lPE13.6)/> 

WRITE (6.9) A.B.C.Zll . Z l 1 P . 2 I 2 . Z l 2P . Z l 3 , Z I 3P . Z I <♦ . Z I b , Z I 6 . Z I 7 
9 FORMAT (//5OX.20HC0MM0N BLOCK /SA VEO / // bX , 7 ( 5X , 1 PE I 3 . 6 ) /SX . 

A 6 ( 5X .1 PE I 3* 6 ) / ) 

WRITE (6.6) ICX.IMX.IPX.IRX.ISX.ITPOS. I TWX , I TX . 1 UX , I XPOS . I YX , iZX. 
t CCX,CMX,CPX,CRX,CSX,CTWX,CTX,CUX,CYx.CZX 

6 FORMAT (//50X.21HC0MM0N BLOCK / LOOK UP / / / l 0 X , 1 2 ( 1 5 . bX ) / 

1 (5X.615X, lpElb.B) ) ) 

IF (IXTAB ,LE. 0) GO TO 100 
IF (IXTAB .LT. 100) 00 TO 22 

13 ■ 95 
60 TO 23 

22 13 ■ IXTAB 

23 1 3 - 10*( 13/10 ♦ t ) 

WRITE (6,7) 

7 format ( //2*fx »77hcommon block /tables/ petab, smtab, tetab, twtab 
1, UETAB, XITAB, YITAB, ZMTAB/ ) 

8 FORMAT (5X, 13. 1P1GE12.S) 

00 2H I » 1,13,10 
X • I ♦ 9 

2R WRITE (6,8) I,(PETAB(J), J « I,K) 

DO 26 1 • 1,13,10 
K ■ 1 ♦ 9 

26 WRITE (6,8) I , ( SM T A B ( J ) , J - I,K) 

00 27 I ■ 1,13,10 

X * l ♦ 9 

27 WRITE (6.8) I,(TETAB(J), J » 1,K) 

00 28 I « 1,13,10 

K ■ I ♦ 9 

28 WRITE (6,8) 1,(TWTAB(J), J a I.K) 

00 29 I a 1,13,10 

X « l ♦ 9 

29 WRITE ( 6 , B ) I.CUETAB(J), J a I,K) 

00 30 I > 1,13,10 

X ■ l ♦ 9 

30 WRITE (6,8) I , ( X I T AB ( 0 ) , J a I,K) 

00 31 I » 1,13,10 

X ■ 1 ♦ 9 

31 WRITE (6,8) 1 , ( Y I T AB ( J ) , J a I , K ) 

00 32 I ■ 1,13,10 

X ■ I ♦ 9 

32 WRITE (6.8) I , ( ZMTAB ( J ) , J ■ I.K> 

WRITE ( 6 , 1 1 ) ICOOL.IOUMP.ITZTAB.AL.COEFCL.CPL.DELXBA.OIaTUB, 

1 FLOWRT ,MaSSL,PRaNOL,RaMOL,RaMOW,REYL.SUM8Ga.SUMOWI , 

2 THICX.TLO.TLl , T L 2 , TLC A « TQL I T E , TUBEN . T W(,C A . ZM VUL , 

3 CPLTAB .RAMTAB .TZTAB.ZMYTAB 

11 FORMAT (//51X.19HC0MM0N BLOCK /C00L///&5X , 1 1 , 3X , I 1 , 3X , I 2/ 

1 11(2x.IP£1Q«R)/11(2X«E10,R>/(10(2X,E11,5))) 

00 33 I a 1,13,10 
X ■ I ♦ 9 

33 WRITE (6,8) l.ULTAB(J), J « I.K) 
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3H 


35 


36 


37 

12 


38 

100 


DO 3H I « 1,13.10 
K » I + 9 

WRITE (6,8) I , (TLTABl J) , J - l .K ) 

DO 35 I « 1,13.10 
K ■ I ♦ 9 

WRITE (6,8) i , (TH1TABC J) * J - I,K) 

DO 36 I ■ 1,13.10 
K ■ l ♦ 9 

WRITE: (6.8) 1, (TWO TAB (J). J ■ I.K) 

DO 37 I * 1,13.10 
K ■ 1 ♦ 9 

WRITE (6,8) l ,(TLCTAB( J) . J * I.K) 

WRITE (6.12) 1 £ X , C E X 

FORMAT ( //SOX , 2 1 HCOMMON BLOCK /NHANCE///3X » I3,6(5X. 1PE15.B) ) 
DO 38 1 » 1,13.10 
K - I ♦ 9 

WRITE (6,8) 1, (ENHTAB(J). J ■ I.K) 

CALL DIRECT 
END 


QUIT HI 
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readin 


c 


c 


c 


c 

c 


c 


c 


c 


2 


HI 1 


SUBROUTINE REAUIN 


REaoQQOI 


COMMON /COOL/ UOOL.IOUMP.ITZTaB.AL.COEFCL.CPL.uELXba.DIaTUb, /COOL/ 

1 FL0*RTiMASSL.PRaN0L.RaM0L,RAM0«.ReTLiSUMQ6A.SUM8I»I , /COOL/ 

2 thick, tlo.tli ,tl2.tlca,tolite. tuber, tagca.zmyul, /cool / 

3 CPL-TAB120) ,RAMTAB( 201 .TZTAB < 20) .ZMYTAB (20) , /COOL/ 

H ALTABI ICO) .TH1TAB! IOC) .TLCTABI 100) ,TLTAB< 100) , /COOL/ 

S T «GTAB (ICO I /COOL/ 

REAL MASSL /COOL/ 


COMMON ZcSeVaL/ NOcTABt IS,R0J,FJG,CJG.0M102,G0GM1 .POMaX.cPO.HO. 
A S0,TCTAb( 20) , CPT AB t 20 ) , BCP ( 20 1 , CCP ( 20 ) , DCP ( 20 ) , 
K STAB <201 .HTAB120) .BARB1 (20) ,BARU2(2Q) ,6 aRB3(20) 


/CSeVAL/ 

/CSEVAL/ 

/CSEVAL/ 


COMMON /input/ ioxmax , I cTAB. IPRINT . ITATAB, IXTAB.MZETA.OXMAX, 
a EPSZ.FJ.G.GAMO.PO.PHU . P I E , PR ANOT , RB AR , ScALE . TO . 

K THETA 1 ,TOLCFA,TOlZET,TOLZME,ZMUO,ZMVIS,ZnSTAN 


/INPUT/ 

/INPUT/ 

/INPUT/ 


COMMON /NHANCE/ 1EX,CEX(6) .ENMTaBI 100) 


/nhance/ 


COMMON /TABLES/ PETABl 1 00 ) , SMTAB ( 100) ,TETAB( 100) »T*T AB { 100) , /TABLES/ 

1 UETAB1 100) ,X I TAB ( ICO) ,T 1 TAB! 100) ,ZMTAB( 100) /TABLES/ 

01 MENS I ON PITABl 100) ,TI TABI100) .TITLE! 13) .VlTABl 100) 
euui valence ipetab ,pi tab i . itetab .titab i , iuetab.v itabi 


namelist /Nam^ 

1 

2 

3 

H 

b 

6 


ALTAB,COEFCL,CPLTAB,CPTaB,OXMAX,ENHTAB,EPSZ,FJ, 
FLOARTtG.GAMO, I COOL, ICTAB, 1 DUMP, IPRINT, I TRTAB, 

itztab* ixtab.massl.hzeTa.po.petab.phu .pi tab. 
r AMOA, r amtab. RBAR, scale, sm Tab, T o, tctab.te TAB, 

THETA I ,THITA8,TITaB,TLTaB,T0LCFa,T0LITE.T0LZET . 

T0LZME.TUBEN.T*TAB.TZTAB,UETA8,VITAB,XITAB,YJTAB. 

ZMTAB.ZMUO.ZMVIS.ZMYTAS.ZNSTAN 


/NaM 1/ 
/nami/ 

/NaMI/ 
/NaM 1/ 
/NaHI/ 
/NAMI/ 
/NAM 1 / 


SCALE ■ 1,0 

READ002A 

MZE1A ■ 7 

REA00027 

ZNS T AN * o.l 

REA00028 

FJ m 770 .2 

REA00029 

G • 32.17H 

RE AD0030 

TOLCFA « l.OE-OH 

tolite « 0,0020 

REA00031 

T0LZME ■ 1.0E-07 

REA00032 

TOLZET - 0.0003 

RE A00033 

OXMAXO ■ OXMAX 

REA0003H 

0 XM A x ■ -28982,0 

READ0035 

REA0(5,2) TITLE 
FORMAT (13A4I 

REA00Q3A 

REALMS, nAMI ) 

IF i 0 XM A X , NE . -28982.0) GO TO H 1 1 
IF (IOXMAX , E9 . 0) GO TO HIS 
OXMAX ■ OXHAXO 

REA0003B 

GO TO H 1 6 

IF (OXMAX ,LE. 0.0) GO TO H 1 H 
IOXMAX ■ 1 

REA000H3 

GO TO H 1 4 

REAOOOHS 
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SIH 10XMAX 
<U5 DXMAX 
<116 


■ 0 

»( ( X I T A B ( IXTAB) - XITAB(l) ) / 100.0 ) • SCALE 

IF (EPS Z .LE . O.O) WRITE (6.7) 

7 FORMAT (//////////56X.I9H*** INFORMATION «*•/ / 30X . HHH 1 . THIS CASE 

i considers two-dimensional flow//30x,3Zh2. the nozzle width is on 

2E FQOT//30X.59H3. THE SIDE WALLS ARE ASSUMED To BE ADIABATIC AND 
31NVISCID//3QX.6BHH. HEAT TRANSFER OCCURS ONLY THROUGH THE ONE F00 
HT WlOE CURVED * ALLS / / 30X , 59HS . THE CALCULATED THRUST LOSS IS BASE 
SD ON TWO CURVED W ALLS// 30X • 6 SH 6 . THE CALCULATED THRUST IS BASED 0 
6 N AN AREA OF ONE BY Z*YR FEE T 7/ 30X . SSH 7 . CHECK THE INPUT VALUES F 
70 R FLOWRT, MASSL. AND TUBEN/I 
WRITE (6.3) TITLE 
3 FORMAT ( 1HI .27X.I3A6//) 

I ERROR ■ 0 

WR I TE ( 6 . 102 ) MZETA 

102 FORMAT! H5H MZETA « VELOCITY PROFILE POWER LAW E XPONENT 27 X 1 H« I H ) 

IF (MZETA .GE. 0) GO TO 2 S 
WRITE (6.300) 

300 FORMAT ( H7H ••ERROR** VALUE MUST BE GREATER THAN ZERO (0). // ) 
l ERROR ■ I 

2G WR1TEI6.103) IPRINT 

103 F ORmAT ( 7 3h IPRINT ■ PRINT AT EVERY CALCULATED POINT(-I) OR AT 
IT INTERVALSCO) « I H I 

IF (IPRINT .EW* 1 .OR. IPRINT • EQ • 0) GO TO SI3 
WRITE (6.S02) 

502 FokMAT ( M5H ••ERROR** VALUE MUST BE ZERO (0) On ONE (I). // ) 

I ERROR > 1 

S13 WRITE16, 10H> IXTAB 
10H F ORMAT ( S2h IXTAB ■ NUMBER OF POINTS 
1 Mb I V ) 

IF (IXTAB .GE. H .AND. IXTAB .LE. 100) 

WRITE ( 6 . 304 ) 

}0H FORMAT (/ 2 X. 10 MH*. ERROR •• VALUE MUST BE GREATER THAN OR EOUaL TO 
1 FOUR ( N ) OR LESS THAN OR EQUAL TO ONE HUNDRED (100). //> 

IEKROR b J 

30 WRlrE(6,10S) ICTAB 

10 S FORMATINSH ICTAB ■ NUMBER OF POINTS IN CP .vs. T TAbLE 27 XlHBl<!) 

IF (ICTAB .EQ. 0) GO TO 37 

IF (ICTAB . GE . J .AND. ICTAB .LE. 20) GO TO 37 
WRITE (6.306) 

306 Format (/2X,9BH»* Error •• VALUE must »t Greater than or equal to 
ITHREE ( 3 ) OR LESS THAN OR EQUAL TO TWENTY (20).//) 

I ERROR ■ l 

3/ WRITU6.I06) ITwTAB 

>06 format ( 73 H itwTau . wall temp, option — aoiabatici--i > 

1*0). TABLE(*l! *I.RI 

IF ( I ABst I TWTAB ) ,EQ. I .OR. ITwTAB .EQ. 0) GO TO S23 
WRITE (6.512) 

512 FORMAT ( 6 7h ••ERROR** VALUE MUST BE ZERO (0). PLUS ONE (1). OR 
1NUS ONE (-1 ) . // ) 

1ERRUR b 1 

523 WRI TE ( 6 , 1 | 1 ) TO 
111 f ORMXT ( *»Bh TO 
1EIS.7) 

IF (TO .GT. 0.0) 

WRITE (6.300) 


READ00H6 

READOORT 


IN X .VS. Y .VS. M 


GO TO 30 


REA00050 
REA0005 1 
REA00052 


REA00055 
RE A 00056 
REA00057 
I NPUREAD0058 
REA00059 


READ0063 

READ006H 


TABLES20X1READ0067 

READ006B 


NEAD0076 

REA00077 

READ0076 


READ0089 
READ0090 
CONSTANT (READ0091 
REA00092 


MIREAU0096 

REA00097 

READ009B 


• ErEE stream stagnation TEMPERATURE 


GO TO Ml 


2 M X 1 H« 


IPREAOOIOI 

REA00I02 
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free stream stagnation pressure 


reaooui 

REAOOl 12 

- STAGNATION RATIO of SPECIFIC HEATS28X|H«1PEi5.READQ1 13 

REAOOl It 

GO TO 47 
GO TO 47 


VALUE MUST BE GREATER THAN One < 1 • 0 > • 


I ERROR * 1 

41 *R1TE(6|112) PO 

112 F ORMAT ( 50H PQ 
1 1 PE 1 S« 7 ) 

IF (PO .gT. O.O) GO TO 43 
WRITE (6.300) 

(ERROR > 1 

43 WRlTE(6,ll3> GAMO 

113 F OHM AT ( 44H GAMO 
17) 

IF (lCTAB .NE. 0) 

IF (GAMO .GT. 1.0) 

WRITE (6.541 ) 

541 FORMAT ( 48H »«ERROR** 

1 ERROR * 1 

47 WRITE (6.115) ZMUO 

115 FORMAT ( 38H ZMUO 
IF (ZMUO .GT. 0.0) 

WRITE (6.300) 

1 ERROR ■ 1 

si write ( 6 a ii 6 ) zmvis 

116 FORMAT (N 7 H zmvis 

115.7) 

WRITE16.1I7) ZNSTAN 

117 FORHAT(nbH ZnSTaN ■ BOUnoaRT LATER interaction 

1.7) 

WR I TE ( 6 . 1 1 8 ) qxmax 

128 FORMAT ( 3 j H DXMAX ■ MAXIMUM STEP SIZE 
IF 1 THE T A l *LT. 0.0) GO TO 44 
WRITE (6.119) THETAI 
1 19 FORMAT ( H9H THETAI 
1 PE 1 5 . 7 ) 


REAOO 1 05 
READO 1 06 
22X1H-READ0107 
RE ADO 1 OB 


// 


STAGNATION 

Go TO 51 


VISCOSITY 


34XIH«1PE1S.7) 


) READO 118 
REAOOl 19 

REAOO 1 26 


» EXPONENT OF V 1 SCOS I T Y -TEMPER ATURE 


REAOO 129 
REAOOl 30 
LaW 25X1H>1 PE READO 131 
REAOOl 32 
READO 1 33 

EXP0NENT27X1H»1PE1SREA00134 
RE ADO 1 35 
RE ADO 1 36 

4 1 X J H» 1 PE l 5 • 7 ) RE ADO 1 37 


INITIAL VALUE OF MOMENTUM THICKNESS 


23X iHi 


1 RE ADO 1 40 
REAOOl 4 I 


WRITE (6.120) PH 1 2 
>20 FORM AT ( r 7 H PHll 
1 IS. 71 

44 WRITE ( 6 , 121) EPSZ 
121 FORMAT ( 5 1 H EPSZ 
1 * 1 PE 1 5 , 7 ) 

IF (EPSZ >CQ. 0.0 
WRITE (6.502) 

ierror » 1 

533 WR1TEI6, 122) rBAk 
122 FORMAT ()X, 3 &HRBAR 
IF (RBAR .GT. 0.0) 
WRITE (6.300) 
IERROR ■ 1 

53 WRITE(6.123) FJ 
123 FORMAT (SlH FJ 


■ INITIAL VALUE OF ENERGY THICKNESS 2 SX 1 H« 1 PERE ADO 1 43 

REAOOl 44 

■ GEOMETRY,., AXISYMMETRIC(»1 . ) , PL AN£ t >0 . ) 2 l X I HRE ADO 1 4 7 

RE ADO 1 48 

•OR. EPSZ .Eq. 1.0) GO TO 533 

RE ADO 1 52 

■ GAS CONSTANT AT STAGNATION, 36X,lHa,lPEl5, 7) 

GO TO 53 

READ01 58 
RE ADO 1 59 

« CONVERSION BETWEEN THERMAL and WURK UNIT521X1HREad0160 


1»|PE15.7) 



IF (FJ 

.GT . 0 

• 0) 


WRITE ( 

6.300) 



IERROR 


* 1 

5b 

WRI T£(6 

.124) 

G 

1 2 H 

FORMAT ( 

57H G 



1 A 1 5X 1 H« 

1 PE 15. 

7) 


IF (b , 

GT. 0. 

0 ) 


GO TO 55 


■ PROPORTIONALITY CONSTANT 
GO TO 420 


IN EOUaTJOn 


RE ADO 1 6 1 


READOl 64 
RE ADO 1 65 
FbM/G*READ0166 
RE ADO 1 67 
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WRITE (A, 300) 

1ERR0R - i 

A 20 WR I T£ ( A , A2 1 ) SCALE 

A2I FORMAT ( 30H SCALE » CONTOUR SCALE FACTOR, A2X, IHw, 1PE1S.7 > 

IF ( TOLCF A ,E0. l.OE-A) 60 TO A02 
WRITE (A.AQt) TOLCF A 

HO I FORMAT { A7H TOLCFA ■ TOLERANCE FOR SKIN FRICTION ITERATION, 25X, 
t 1 H» , 1PE1S.7 ) 

AQ2 IF (TOLZET ,EQ, 0.0003) 60 TO AOS 

WRITE (A.AOA) TOLZET 

ROR FORMAT ( R9H TOLZET ■ TOLERANCE FOR SHAPE PARAMETER ITERATION. 

1 ?3X, 1H-, l PE 1 5. 7 ) 

ROS IF (T0L2ME ,E0. l.OE-7) 60 TO 205 

WRITE ( A • R07 ) TOLZME 
RQ7 FORMAT ( ASH TOlZME ■ TOLERANCE FOR MACH NO. - TEMPERATURE REL AT I ORE ADO 1 83 
IN ITERATION, 7 X , |H«, IPE1S.7 ) REA0018R 

205 WRITE (A. 900) I TZT AB , J DUMP .FlOWRT , MA SSL , R AMD w .C oEFCl , TUbEN . 

i tol i te , i cool 

VOO FORMAT (IXiAShITZTAB • NUMBER OF POINTS IN T .VS. CPL .VS. RAM0L . 

IVS« ZMYUL TABLES . 3X , IH* . 1 R/ 1 X .63HI0UMP ■ COOLANT FLOW OPTION — S 
’A ME 0JRECTI0N(»1 ) , REV£RS£I»0),8X,1H«,1r/1X # S2HFL0*RT ■ COMBUSTION 
3 CHAMBER MASS Flow RATE ILBM/SEC) . 19X, 1H«. 1PE1S.7/1X,AIhMASSL ■ C 
roolant mass flow rate (lbm/sec) ,3ox,ih*.ei5.7/ix.rahramdw » heat 

SCONOUCT JVITY OF THE CHAMBER WALL ,2SX , 1 H» , E 1 B , 7/ 1 X .31HC0EFCL ■ COEF 
oF I C I ENT OF COOLINOiROX. i H- ,E 1 S . 7 / 1 X . 20HTUBEN ■ TUBE NUMBER. 5lX, 

7 1H»,EIS.7/1X,52hT 0LITE » TOLERANCE FOR TOTAL HEAT TRANSFER ITERAT 
8t0N,!?X,lH«.El5.7/lX,ARHIC00L » COOLINS OPTION — WITH C00LIN6(«l 
9), WITHOUT CUOLING(-O) ,7X, IH-, IR I 
IF IlCTAU ,LE. 0 .AND. 1 TZT AB ,LE. 0) 60 TO 11 

WRITE (A, 131) 

131 FORMAT (//2X.1MI , 5 X , l 3HSPEC I F I C HE AT , 5X , H H TEMPER ATURE . SX , 

1 1 2HC00LANT TEMP ,SX , IOHCOOLAnT CP.5X,12HC0NDUCT1VITY,SX. 

2 9HV1SC0SITY) 

IMAX • aMAXKICTAB.ITZTaB) 

DO 133 1 - l.IMAX 

IF (I ,LE. 1CTAB .ANU, 1 ,LE. ITZTAB) 60 TO 1 3o 
IF IICTab «6T. ITZTAB) 60 TO 1 32 

WRITE (6.1) I .TZTABt I ) .CPLTABI I > ,RAMTaB< I ) .ZMYTaBI I) 

1 FORMAT (I3,R1X,F9,3,AX,F10,A,SX,F12.10i3X.F12.1O) 

60 TO 133 

130 WRITE ( 6 . R ) I ,CPTAB( I ) .TCTABI l » ,TZTAB( l ) .CPLTABU ) .RAMT aBI I ) . 

1 ZMYTAB(I) 

R FORMAT (13,SX,Fl3,lC,6X,F9.3,eX,F9.3,6X,F10.6,5A,Fl2.l0.3X,Ft2>10) 

60 TO 133 

132 WRITE (A.S) I ,CPTAB( I ) , TCTABI I ) 

5 FORMAT ( 13.SX.F13. 10.AX.F9.3) 

133 CONTINUE 

IF UCTaB «LE. 0) 60 TO It 

II ■ ICTAB - I RE ADO l 93 

DO 59 I • 1 , I I REA0019R 

IF ( TC T aB ( 1 + 1 ) .6T, TCTAB(I)) 60 TO S9 

WRITE (A. 310) 

310 FORMAT ( / 2X » 99H* • ERROR •* TABLE OF SPECIFIC HEATS - TEMPERATURE V 
1 ALU£S must be IN MONaTOnICALLY I NCREAS IN6 ORDER.//) 

1 ERROR . 1 READ0201 

s 9 CONTINUE REA00202 


RE ADO 1 70 
REAOO 171 
REA00172 


REAOO 175 
RE ADO 1 7 A 


REAOO 1 79 
REAOO 1 BO 
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READ0205 


.IF (TCIABIU j«T. 0.0) 60 TO 61 


«R1TE(6,3|2) 

31? FORMAT (/2X.87H** ERROR •• TaBLE OF SPECIFIC HEaTS - TEMPERATURE V 
1ALUES MUST BE GREATER THAN ZERO (0).//) 


I ERROR » 1 

61 00 63 I a I.ICTAB 

IF (CPTABII) .GT. 0.0) 
WRITE(6,3l3) 

313 FORMAT l/?X.8?H*» ERROR 
I VALUES MUST BE GREATER 
1ERR0R • 1 

63 CONTINUE 
11 00 6S I * I.IXTaB 

if izmtab(I) .gt. o.o) 

WRITE (6.319) 

31*1 FORMAT (/2X.97M** ERROR 


GO TO 63 

** table of specific heats 

THAN ZERO (0).//) 


REA0021 1 

SPECIFIC heat 


REA00219 

READ021S 


GO TO 65 

•* table of mach number distribution - mac 


1H NUMBER VALUES MUST BE GREATER THAN ZERO (0).//> 
IERROR » 1 

65 CONTINUE 

II - IXTAB - 1 

00 67 I - 1, II 

IF (XITAB(I*1) .GE. XITAB(D) GO TO 67 


READ0223 

REA0022H 

REA0022S 

READQ226 


WRITE (6.316) 

316 FORMAT ( 90H ••ERROR** TABLE OF CONTOUR DESCRIPTION. / 69h AXI*L DREA00229 

,ist a nce values t x i must be in monotonically increasing order. // ireaoo23o 
IERROR ■ 1 REA00231 

67 CONTINUE READ0232 

IF(ITWTaB) 19,13,12 REA00233 

. 12 Oo 69 l • 1, IXTAB REA0Q239 

IF (TWTaB(I) .GT. 0.0) GO TO 69 
WRITE (6.317) 

317 FORMAT (/2X.102h*« ERROR •• TABLE OF WALL TEMPERATURE DISTRIBUTION 

i - temperature values must be greater than zero ioi.//) 


69 

IERROR • 1 

CONTINUE 


READ023? 

READ0240 

13 

GO TO 19 

IF (TWTABd) .GT. 0.0) 

GO TO 19 

KEA0024 1 

14 

WRITE (6.317) 

IERROR b 1 

IF (SCALE • EG • 1.0) GO 

TO 929 

READQ24& 


00 923 I b I.IXTAB 
XITABI I | b X I T AB ( I ) 

• scale 

READ0250 

423 

Y I T AS ( I ) » Y l T AB ( I I 

• SCALE 

REA00251 

424 

IF (1TWTABI 137,135. 190 




135 write (6.136) TlATAB(l) 

136 format ( //2X . 1 BHWALL TEMPERATURE b.FZO.B) 


137 WRITE (6,6) 

6 format (ihd 


WRITE (6,138) 

138 FORMAT (3X.1H1»1ZX.SHaXIAL.11X,6HRaDIAL.10X,HHMaCH.9X,8hPHES5URE. 

1 9X, l 1HSTATIC TEMP.7X.8HVEL0CI TY,6X,9 hM0LECULAR/ISX,6H(FEET) . 1 IX, 

2 6H|FEET),9X,6HNUrtBER.8X,8H(LB/FT2),9X,HH(0EGREE5 R),7X, 


3 8H(FT/SEC) .8X.6HWEIGHT ) 

WRITE (6.139) (1,X1TAB(I).YITAB(1I.ZMTAB(I),PETAB(I).TEtAB(I). 

I UETABI I ) ,SMTAb( I ) , I ■ I.IXTAB) 

139 FORMAT <I9,6X,Fll.6,6x.F>l.6,6X,F9.6,6X,Fl0.3,6x,F9,3,6x.F9.3.6x, 
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1 F9.6) 

GO TO lib READ0260 

1 HO WRITE (6.6) 

WRITE (6.112) 

112 FORMAT <3X , 1HI , 12X .SHaXI A l, l 1X.6HRAOI AL, 10X,1 HMaCH,9X,8hPRESSURE, 

1 1X.I1HSTATIC TEMP.7X,8hVELOCITY,6X»9hM0LECOLAR,6X,VH«ALL TEMP/ 

2 1BX.6H(FEET).IIX.6H(FEET).9x.6hNUMBER.8X.Bh(LB/FT2),1X, 

3 1 iMtOEGKEES R ) . 7 X t 8H ( F T / SE C ) . 8 X , 6M WE 1 GH T . SX , l I H ( DEGRE ES Rll 

WRITE (6.111) ( 1 .XlTAB ( I ) ,Y1 TAB! I ) ,2MTAB( I ) .PETA8 t I ) ,TETAB( I ) . 

I UETAB ( I ) .SMTAB ( l ) .TWTA8 ( I ) . 1 ■ l.IXTAB) 

11 1 Format (11.6X.F11.6,6x»F11.6,6X.F9.6.6X,F10.3,6x.F9.3.6x.F9.3.6X. 


1 F9.6.6X.FV.3) 

115 IF (ICOOL .GT. 0) WRITE (6.33) ( I , AL T A8 ( I) . T L T A B (I) . T H I T A B ( 1 ) » 

1 ENHTAb ( I) , I ■ l.IXTAB) 

33 FORMAT ( IHJ ,bX , 1HI ,&X , 1 7HC00EANT TUBE A RE A , 5 X , 1 VHCOOL A N T TEMPERATU 
1RE.5X.11HWALL THICKNESS, 5X, 11 HE NHANCEM£NT/16 x,13h(SuUaRE FEET) i 7 X , 

2 17m(DEGREES RAnKINE) . 13X.6H<FEET1 .7X.7HFACTORS/C 1X.13. 1 1X.F1 1*8, 

2 1 IX .F 1 0.3.8X ,F 1 1 .8 ,5X .F 1 l .8) ) 

IF (IXTAB >LE. 1) GO TO 260 
00 257 I ■ 2. IXTAB 

IF IXITAB(I) .GT. XITaB(I-I)) GO to 257 
WRITE (6,212) 

212 FORMAT I 33H •*ERROR»» TABLE OF XlTAB VALUES. // ) READ0269 

1ERR0R ■ 1 READ0270 

257 CONTINUE REA0Q271 

260 IF (ITWTAB .LT. 0) GO TO 77 

IF l THE TA I .GE* 0.0) GO TO 77 
WR1TEI6.76) 

76 FORMAT ( // 99H *»ERROR»« MACH ONE START DOES NOT PRODUCE REAS0NABREAO0276 


1 LE VALUES FOR OTHER ThAN AN ADIABATIC WALL CASE. 
1 ERROR ■ 1 

7/ IF (IERRUR .LE. 0) RETURN 
CALL QUITS 
END 


// ) 


READ0277 

READ0278 


READ0281 


94 


SEVAL 


SUBROUTINE SEVAL ( INOl .AA.BB.CC) 

c 

COMMON /CSEVAL/ N0cTAB,IS,KOJ 1 FJG,CJG l 6Ml02 t OOSM|,POMAX.CP0,HO. 
a S0.TCTAS120) ,CPTAB(20I .8CP(20) ,CCP<20) ,0CP(20> i 

K GTAB(2Q).HTAB(20).BARBI(20).BARb2(20),BARB3(20> 

C 

COMMON /INPUT/ 1 DXMAX , l CTAB , I PR l NT , I TAT AB , I XT AB ,M/ET A .DXMAX , 

A EPSZ.FJ.GiGAMOiPO. PH II. PIE. PRAN0T.RBAR.SC ALE. TO t 

K THETA I .TOLCFA . TOl ZET , T OL ZflE , ZMUO , Z M V l S , ZnS T AN 

C 

C DEFINE The FUNCTION ROUTINE TO BE USED BY SEVAL 

C 

GaPFIT.G.AI .81 .Cl .Oil ■ G*A 1 *ALOG ( T ) ♦ ( 8 1 ♦ ( C 1+0 1/3. »T I «T) »T-PR 
T * AA 
A > BB 
B » CC 

IF ( lNOI-2) 3. 1 . 1SS 
I B « B/FJG 

IF (1CTAB .GT. 0) GO TO 3b 
T ■ B/CPO 
A « CPO 
GO TO 600 

3 IF I INOl .LT. 1) GO TO 10 
IF (ICTAB .GT. 0) GO TO 10 
A = CPO 
B * C JG* T 
GO TO 600 

Ibb PR * ROJ*ALOG( A/POMAX) 

160 STAB * GAPP(TCTAB( ISI .GTABI lSl .BARB1 ( IS) ,BARB2( |S) .BARB3I IS) . 

A OCP(IS)) 

IF IB .GE. STAB) GO TO 17B 
IS * IS - 1 

IF (IS .LE. 0) 60 TO 17 

GO TO 160 

17b STAB * GAPF(TcTAb( IS»1 ) .GTABI IS) .BARBl US) .BARB2I IS) .BARB3I IS) . 
A DCP(IS)) 

IF (B ,LT. STAB) GO TO IVO 
IS « IS * 1 

IF (IS .GE. 20) 00 TO 17 

00 TO I 7b 

1V0 IF (IS .GE. noctab .or. is .LE. O) GO TO 17 
T TP - TCTAB(IS) 

FP « GAPF (TCTAB( IS ) ,GTAB( IS ) .BARBl ( IS) ,BARB2( IS » ,BARB3( IS ) , 

0 DCPIIS)) 

TTPP ■ TLTABI IS*1 ) 

FPP • STAB 
GO TO 75 

10 IF IT .GE. TCTAB(IS)) GO TO 15 
IS- IS- 1 

IF (IS .LE. 0) 60 TO 17 

GO TO 10 

15 IF (T .LT. TCTABIIS+DI GO TO 16 
IS ■ IS + I 

IF l IS .GE. 20) GO TO 17 
60 TO lb 


/CSEVAL/ 

/CSEVAL/ 

/CSEVAL/ 

/INPUT/ 

/INPUT/ 

/INPUT/ 


SEVA 16 
SEVA IV 
SEVA 20 
SEVA 21 
SEVA 22 


SEVA 27 


SEVA 32 


SEVA 38 


SEVA S3 


SEVA 50 
SEVA b 1 


SEVA 5S 


SEVA 57 


95 



SEVA 59 


16 IF (IS .GT. 01 GO TO 19 
17 6RITEI6, 18) IS.INDI .T.A.B 

18 FORMAT (I7H0SEVAL F A It ORE , , • . 5* . 9H I S a , I 9 . 5X , 6H I NO 1 «.I2.SX, 

A 3HT *. 1PEIS.7.5X.3HA ■ . 1 PE 1 5 . 7 • 5X . 3HB a.lPElS.7) 

CAUL QUITS 

19 IF (IS . GE . NOCTAB) GO TO 17 
IF ( 1 NO 1 ) 70.65,60 

60 OELT at- TCTAB(IS) 

B - HTAB(IS) ♦ CPTABllSjaDEUT ♦ 0 . S*BCP ( I S ) *OELT »»2 ♦ 

A CCP(IS)/3.0»0ELT»O ♦ 0.25*0CP( IS)*DELT**8 
B « B*F JG 
GO TO 191 

65 PR ■ R0J»AL0G(A/P0MAX) 

B ■ GAPF(T,6TAB( IS) .BARBl I IS) ,BARB2( IS) .BARB3I IS) ,OCP( IS) > 

GO TO 600 

70 A » POMaX*EXP( ( GTaB ( IS) ♦ B A RB 1 ( I S ) • AL Ofi ( T » a (BARB2IISI a 
A ( BARB3 (IS) a DCP( 1 S I / 3 . 0*T ) • T ) *T - B)/ROd) 

GO TO 600 

35 IF (B .GE. HTAB(IS)) GO TO 50 
IS ■ IS - 1 

IF (IS ,LE. 0) GO TO 17 
GO TO 3S 

50 IF (B .LT. HT AB ( I S+ 1 ) > GO TO 51 
IS « IS a l 

IF (IS .GE, 20) GO TO 17 
GO TO 50 

51 IF (IS , GE , NOCTAB .OR. IS <LE« 0) GO TO 17 
TTP • TCTABI IS) 

FP a HTAB(IS) 

TTPP ■ TCTABI IS+1 > 

FPP » HTAB ( IS* 1 ) 

7 S T TO ■ (TTP • (FPP - B ) - TTPP • (Fp - 8 )) /( FPP - FP) 

IF (INDl .GT. 2) GO TO 21S 
OELT a TTO - TCTABI IS) 

FO ■ HT AB < IS) ♦ CPTABl ISIaDELT a 0 . 5*BC P ( 1 S ) *0 ElT*« 2 a 
A CCPI IS)/3.0*0ELTa*3 ♦ 0 . 25 aOCP ( I S ) aDELTa *H 

GO TO 220 

215 FO ■ GAPF ( TTO.GTABl I S) .BARB1 ( IS) iBARB2< IS) ,BARB3( IS) .OCP( 1 S) ) 
220 TTlP ■ (TTO • (FPP - B ) - TTPP *(F0 “ B ) ) /( FPP - FO) 

TT1PP a (TTO a ( FP - B ) - TTP a ( FO - B))/( FP - FO) 

N ■ -1 

TAU a TTO 

SF a FO 

109 IF (ABSItSF - B>/B> « LE . 1.0E-7) GO TO 100 

IF (SF ,LE. B) GO TO 135 
TTPP a TAU 
FPP a SF 
GO TO 130 
100 T a Tau 
GO TO 225 
135 TTP ■ TAU 
FP a sF 

130 IF ( N )115 , 120 , 125 

1 1 5 N ■ 0 

TAu • TTlP 
GO TO 65 


SEVA 61 


SEVA 68 


SEVA 71 


SEVA 79 


SEVA 77 


SEVA 80 


SEVA 87 


SEVA 92 

SEVA 99 
SEVA 95 

SEVA 97 
SEVA 98 


SEVA 101 
SEVA 102 
SEVA 106 
SEVA 107 
SEVA 106 
SEVA 109 
SEVA 110 
SEVA 111 


96 


SEVA 112 
SEVA 113 


120 N ■ l . 

TAU ■ TT1PP 

6S IF (TAU .LE. TIP .08. TaU .GE. TTPP) GO TO 130 
IF ( 1 NO 1 .LE. 2) GO TO ?S 

SF - G*PF t TAU.GTABUS) «BARBM1S> .BA8B2.US) *BARB3C IS) .DCpl IS1 > 

GO TO IOH SEVA 118 

95 DELT ■ TAU - TCT A8 1 IS )- 

SF » HTAB(JS) ♦ CPTABl IS)*DElT ♦ 0 . 5*BCp ( I S ) *DE|_ T • • 2 ♦ 

A CCP(IS)/3.0«0ELI«»3 * 0.2S*DCP ( I S I *DELT • *H 

GO TO 10** SEVA 122 

1 25 IF UtFPP - FPl/tFP ♦ FPPJ) .GT. 0.0010) GO TO 75 
T « l T TP* ( FPP - 8) - TT'PPMFP - B I > / ( FPP - FP) 

225 IF IIND1 .'GT. 2) GO TO 600 
DEL T « T - TCTABIIS) 

141 A ■ CPTabIIS) ♦ BCP(1S)*DELT ♦ CCP(IS)*DELT**2 ♦ DCPl IS)*0ELT**3 
600 AA * T 

BB « A SEVA 130 

IF (1ND1 .EQ. 2) RETURN 
CC ■ B 
RETURN 

ENO SEVA 134 


I 
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START 


subroutine start 

COMMON /CO FllF/ I F I N T , AF I N T . 6 F I N T , C F I N T , MM I N T , Tf I N T 

COMMON /CSEVaL/ NOCTABilS,ROj,FJG#CJG*OMl02t6OGMl|POMAX % cP0tHO» 
S0iTCTAB(20) ,CPTAB(20) » BC P ( 2Q ) ■ C CP ( 20 ) » D C P ( 20 ) * 
G T AB ( 20 ) * H T A B t 20 ) ,8ARBl(20) »BAR&2(20) * B a R B 3 ( 20 ) 


COMMON /INPUT/ 


A 

K 


DXMAX . 
SCALE ,TOt 
2NSTAN 


common /INTER/ 


A 

K 


COMMON /LOOKUP/ 


1 

2 


IDXMAX.ICTAB.iPRlNT.ITwTAB, IXTAB. MZETA, 
EPSZ.FJ.G.GAMO.PO.PHI I , P I E • PR AND T , RB aR , 

THETA I fTOLCFA. TOlZET. TOL Z ME , ZmUO . Z M V I S , 

CFAGT ,CFAGP tCHPARl # DX .OXRHO.HE ,HA # IBEG.M 2 ETAM # 
OOMZET .PhIP,PRE103#RH0E«RH0UE,RMZETA , THE TAP • 
XIBASE.XIEND .ZETATM .ZMZETA .ZMZETh, ZMZETP 

ICX.lMX.lPX.lRX.tSX.lTPOS.lTKX.iTXiiUX.lXPOS.lYX 
1 Z X * CC X I 6 ) , C M X ( 6 ) * CP X ( 6 ) . CN X ( 6 ) • C SX ( 6 ) » C T W X ( 6 ) i 
C T X ( 6 ) . C UX ( 6 ) » C Y X ( 6 ) » C ZX ( 6 ) 


'COMMON /OUTPUT/ 


BDELTA.CF,CH*DELTA,OELSOT,t)ELSTR, FLAT .FORCE. HG» 

PE. PH I .QW.SUMQDa.TE.THETA . T ft . UE . X . XLA RC .YR.Zl , 
Z2>Z3.ZM»Zb.ZETA.ZME 

/SAVED/ A.B»C.2I1.ZI1P.ZI2.ZI2PiZI3,ZI3P,ZIH.ZIS,ZI6,ZI7 /SAVED/ 

/TABLES/ PETAB<t0Q>.SMTAB<100>.TETABU00>.T#TABUQ0). 

UETAB( 100) .XITABUCO) .Tl TABU 00) .ZMTAB( 100) 


STAR 1 

/C Of I IF/ 

/CSEVAL/ 

/CSEVAL/ 

/CSEVAL/ 

/INPUT/ 

/INPUT/ 

/INPUT/ 

/INTER/ 

/INTER/ 

/INTER/ 

/LOOKUP/ 

/LOOKUP/ 

/LOOKUP/ 

/OUTPUT/ 

/OUTPUT/ 

/OUTPUT/ 


COMMON 

COMMON 


I 


I M X « -1 

ITS - - \ 

IT*X « -1 

l m 0 

SI ■ I ♦ 

IF { 2 M T A B ( I ) - I. 

10 IF (I .LT. IXTAB) 
25 WRITE (6.1000) 

1000 FORMAT </3&X.62h** 
1LU0E M » I • 0/ / > 


/tables/ 

/tables/ 


STAR 35 


1 

) 10 . 

GO TO 


20 . 

S 


IS 


star 

STAR 

STAR 

star 


37 

38 

39 
HO 


start failure ** mach number table does not inc 



CALL 

quits 





20 

X 

m 

X I TAB l I ) 


STAR 

M3 


1 BE Q 

u 

1 + 1 


STAR 

MM 


60 TU 

50 



STAR 

MS 

IS 

IF (l 

• LE • 1) 

GO TO 2S 





XG « 

X 1 TAB ( I ) 






ZME 


ZMTAB ( l ) 


star 

M 9 


X 


* 1 TAB C l - 1 ) ♦ 

( X l T A 8 ( l ) - XITAB( 1-U 1 

/ ( ZMTAB ( I )STAR 

SO 


1 


- Z M T A B ( I - 1 ) 

) • ( 1 • - Z M T A B ( 1-1 ) ) 

star 

51 


J 


0 


star 

52 


IBEG 


I 


star 

53 

35 

J 


J + 1 


star 

5 M 


XQ 


XG 


star 

55 


ZMU 


ZME 


star 

56 


XG 


X 


star 

57 


CALL 

xnterp ( 

X. ZMEi ZmEp, 

I M X . X I T A 0 * ZMTAB, iXTAB, 

Cmx, imx ) star 

5B 


98 



IF (ABScZME - 1*0) »U£t T0LCFA) SO TO SO 
ZMX « US - XOI/UME - ZMOI 



X 

m XO ♦ II. • ZMO 1 * ZMX 


STAR 

61 


IF (J .LE. SOI GO TO 35 





• RITE (8.1010) 




1010 

FORMAT { 6 R H 

••start failure. . . mach no. calculation eXceeoeo So 

ISTar 

25 

iterations. // l 


star 

26 

50 

call xmterp 

e x, ZME. ZMEP, IMX, xitab. ZMTAB, ixtab. cmx. 

IMX ) 

star 

69 


CALL XNTERP 

U.TR.YRP.IYS.XITAB.YITAB.IXTAB.CYX.IMX) 





CALL GETPT (ZKE t PSE f TE) 





CALL SEVAL ( 

l. TE. CPE. HE > 


star 

67 


OA ME 

* CPE / | CPE - RBAR / Fj 1 


star 

68 


He 

• HO - HE 


star 

69 


UE 

a SORT! 2.* HB ) 


star 

70 


RHSE * PSE/TE/KBAR 





ZKO 

■ ZMUO * ( TE / TO ) «* ZMVI5 


star 

72 


haw 

* HE +CPRANOT ** < !•/ 3. )> ♦ HB 


star 

71 


CALL SEVAL C 

2. TAW, C P A W * HAW I 


star 

7 9 


IF l ITWTAB 

) SS, 6C, 6S 


star 

75 

se 

HA 

* HAW 


STAR 

76 


TA 

* TAW 


star 

77 


60 TO 70 



star 

78 


K« a TWTAB ( 21 





TW 

« TWTABi 1 I 


star 

80 


GO TO 70 



star 

81 

65 

CALL XMTERP 

{ X, TW, TWP, I T W X , XITAB, TWTAB, IxTAB, CTwX, 

IMX ) 

star 

82 


call seval ( 

It TW, CPW, HW 1 


star 

83 

70 

AF tNT 

* HW 


star 

89 


BP I NT 

■ HO - HW 


star 

85 


CF I N-T 

* - HB 


star 

86 


TF INT 

« TE 


star 

87 


mmint 

» HZETA 


star 

86 


IF I NT 

= 1 


star 

89 


CALL II9TZET 

C 0* • 1 . , 21 1 ) 


star 

90 


IF INF 

a 2 


star 

91 


CALL jmtzet 

( Q. . I.. Z12 ) 


star 

92 


t>ELSOT 

« ( I* / ZHZE T A - Z12 1 ■ / III 


star 

91 


IF fEPSZ *EQ 

• 0 • Q 1 Go TO 72 





ERASES * YRP/TR 





*»0 TO 71 



star 

96 

72 

ERASES 

X 

o 


star 

97 

71 

ERASE* 

« < 1. + DELSOT ) / 1 i. ♦ 1 GAME - 1 * I / 2. 

> * 

star 

98 


l 

zme p ♦ erases 


star 

99 


IF < ERASE * . 

NE. 0.0 GO TO aO 





• RITE (8.1320) 

1020 format e/zx.7NH«* start failure *• initial values for phii and the 

lT*l CANNOT BE C0HPUTED./3X ,S5H* CHECK SLOPES OF MaCH NUMBER AN» CO 


2NT0UR INPUT TABLES.//! 

BO ERASE! « 17.2 • < TO “ TA* > / T*M STAR IQ3 

ERASE2 « 30s . • ( TE - TO 1 / TA* STAR 10B 

CTHET * C. 5C»S8RT ( t .0 ♦ TRP*»2 I /ERaSEB 
CRT2 * ( TAA/TE I *•( 1 .0 - ZM V I S I *RHSE*UE/ ! ZMU*C THET I 

ERASE3 • TAA / TE STAR IDS 

CFA a .001 STAR 109 

JE mQ STAR HO 


85 JE 


m JE ♦ 1 


STAR 111 




CFG > CF A 

star 

112 


THETA * CFG • CTHET 

STAR 

113 


CR * CRT2 • ThETA • THETA 

STAR 

11H 


CFB « CFEVAL(CR) 

TTA« ■ i.o ♦ ERASEI*SQRT(CFB/2.0) ♦ ERASE2»CF8/2.Q 
IF (TTA* .GT. 0.0) GO TO 120 
C F A « 3 . 0*CFG 

STAR 

US 

105 

ZH m CF A 

STAR 

11? 

1 10 

Z2 > CFG 

IF (JE • LE . 50) GO TO 85 

STAR 

120 


WR I TE ( 6 . 1030) 

STAR 

123 

1030 

FORMAT (/2X.7HH.. START FAILURE •• INITIAL VALUES FOR PhII ANO 
1 T A 1 CANNOT be C0MPUTED»/3Xi27H* CHECK SKIN FRICTION OaTa.//) 

THE 



60 TO mo 

STAR 

12H 

120 

CF A » CFB / ( ERASE3 • TTA* •• ZMVIS ) 

IF (ABSUCFA - CFGI/ICFA ♦ CFG)) ,LE. TOLCFA) GO TO IHO 
IF (JE .LT. 2) GO TO 105 
23 « ZH 

STAR 

125 


21 >22 

STAR 

129 


2** > CF A 

STAR 

130 


22 > CFG 

2S5 » ( ZH - Z3 ) / ( Z2 - Z 1 ) 

CFA > (ZH - ZSS*Z2 > / ( 1 .0 - ZS5) 

STAR 

131 

1 HO 

GO TO no 

THE T A I • CFA*CTHET 
PH 1 I > THETAI 

STAR 

13H 


CFAGT > CFA 

STAR 

138 


2ETA > 1. 

WRITE ( 8 1 1 OHO ) X, YK, THETAI , PHI I 

STAR 

139 

iqho 

FORMAT ( 9HH0INIT1AL VALUES FOR ENERGY ( PHII) AND MOMENTUM < 

ThETSTaR 

30 


1 A I ) THICKNESSES CALCULATED AT THROAT... / SH X a I 1PE1H.7. 5 X, STAR 

31 


2 HHY > . E1H.7. SX. GhThETAI >. ElH.7. 5X. GhPHII ■ . E|<*.7 // 

> STAR 

32 


RETURN 

STAR 

1 H 1 


END 

STAR 

1H2 
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. ANTERP 

SUBROUTINE XnTERP (X.Y.YP.IXIN.XAR.YAR.IAR, CAR, IPOS) 

C 

DIMENSION C 16) .CARC6) ,XAR| IAR1 . X I ( 4 1 , T AR ( I AR I . T I ( 4 1 
C 

IXO-IXIN 
I XMAX»I AR-l 
1X-IPOS 

PO li I ■ U6. 

11 c ( I ) “C ar i 1 1 

IF ( 1X0 *6t* 0) NO TO 13 

12 IFIRST*) 

IXO*1XHaX* 2 
1 X«l 

13 IF (IX «I.E« 01 60 TQ 12 

20 IF (X ,6E. XAR(IX)) 60 TO 25 
IX ■ IX » l 

IF (IX »GT* 01 60 TO 20 

22 *RlTE(6»23> X. XAUU). XAR( iXMAXt 1 ) , YAR(l), YAr ( I XMAX* I ) 

23 FORMAT ( 2BH0XNTERP OUT OF RAN6E.., X *, IPEI5.7. 8H , XII) 

1 El5*?t BHi XtN) «. EIS.7 / 43X, 8h. Y(l) *. E1S.7. 

2 BH Y»N> EIS.7 // ) 

CAUL QUITS 

2S IF (X .UE. XAR ( I X* 1 1 1 • GO TO 27 
IX * IX ♦ 1 
I F I I X»I XHAX ) 2S.2S.22 

27 PO 28 1*1.4 
l 1 “ l X-2* I 

xim-XARim 

28 Y I I I ) *Y AR ( I 1 ) 

0X2 * X 9 XI (2 1 

0X32-XI fil-XI (2) 

IF (IX - 1X01 40*31*60 . 

31 1X060-0 

IF (IX .61. 1 > 60 TO 33 

32 l 60s- I 

60 TO 101 

33 IF (IX .(Li. I XHAX) 60 TO 3S 
IF ((FIRST *C«t 01 40 TO >1 

(FIRST • O 
(60*1 
6q TO 45 

34 I 60» 1 

60 TO 100 

35 160*0 

60 TO 100 
40 1X060*91 

IF (IX * LT* IXO - I) 60 TO 42 

c 1 4 1 - cm 
6(51 *C (21 
C(4|*C(3> 

60 TO 43 
42 Cr4)*YIl2I 

DX42*XI 141-XI (21 
PY32«Yl(31-rl(2) 

D VOX 32*0 V 32/OX 32 


XNTE 4 
XNTE S 
XNTE 6 

XNTE 21 

XNTE 23 
XNTE 24 
XNTE 25 


XNTE 30 
XNTE 31 
XNTE 32 
XNTE 33 
XNTE 34 


XNTE 37 
XNTE 38 
XNTE 39 
XNTE 40 
XNTE 41 

XNTE 49 

XNTE SI 

XNTE S3 
XNTE 54 


XNTE 58 
XNTE 59 
XNTE 60 
XNTE 61 
XNTE 62 
XNTE 63 
XNTE 64 


XNTE 67 
XNTE 68 
XNTE 69 
XNTE 70 
XNTE 71 
XNTE 72 
XNTE 73 
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C(6)a(0yoX32-( Yl (Ml-YI (21 1/DXM21/I XI (3>-Xl <*»> » 
C(5)»DY0X32-C(6)»DX32 
IF (1XOGO .GT. C) GO TO 100 
H 3 IF (IX .LE. 11 GO TO 32 

IGO a Q 

MS C( 1 » -V I ( 1) 

0X2 1 I (21-XI < 1 ) 

0 X 3 l »X 1 ( 3 I “X l ( 1 ) 

0Y2 1 «Y I (21-YI ( 1 1 
OYOX21»DY21/UX2l 

C(3>»(OYOX2l-(Yi (31 -Y1 1 1 ) 1/0*31) /.I XI (21-XI <3! 1 
C(2)»DY0X21-C(3)»0X2I 
IF(IXOGO) 100,100,42 
60 IXOGOsl 

IF (IX .GT. 1X0 ♦ 1J GO TO MS 

C( 1 1 » C ( M 1 

C(2)>C(S) 

C(3)»C(6I 

62 IF (IX . GE • IXMAX1 GO TO 34 
IGO a 0 
GO TO M2 

100 0X1 > X - XI ( 1 ) 
rei»(C( 3 i«oxi*c (2i)*o xi ♦cm 
YP8i*C(31/»5*0Xl*C(21 

IF (IGO .GT. 0) GO TO 110 

101 YB2 a ( C ( 6 ) *0 x2 ♦ C(5))*DX2 ♦ C(M) 
YPB2*C(4)/.5»0X2*C(SI 

IF (IGO .UT. 0> GO TO 120 

U1 ■ DX2/DX32 

U2*0 1 *U 1 

U3*U2»U 1 

A1*3,»U2-2*»U3 

AlP*6.»(Ul-U2)/0X32 

Y ■ ( l ,-A 1 )*YBI*A1*YB2 

YP*(l.-Al)*YPBl-AlP*(Y8l-YB2)*A|*YPB2 
105 I X I N* I X 

IF (IXOGO • EQ « 0) RETURN 
00 107 I ■ 1,6 
107 C AR ( I ) »C ( I ) 

return 

110 Y*YB1 

Y P* Y PB 1 
GO TO 105 
120 Y ■ Y B2 

YP-YPB2 
GO TO 10S 
ENO 


XNTE 

79 

XNTE 

75 


XNTE 

7? 

XNTE 

60 

XNTE 

81 

XNTE 

82 

XNTE 

83 

XNTE 

8 9 

XNTE 

85 

XNTE 

86 

XNTE 

87 


XNTE 

90 

XNTE 

91 


XNTE 

9** 

XNTE 

97 

XNTE 

98 


XNTE 

102 

XNTE 

106 

XNTE 

107 

XNTE 

108 

XNTE 

109 

XNTE 

110 

XNTE 

111 

XNTE 

112 


XNTE 

117 

XNTE 

119 

XNTE 

120 

XNTE 

121 

XNTE 

122 

XNTE 

123 

XNTE 

129 

XNTE 

125 
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zetait 

SUBROUTINE ZETA1T 
C 

COMMON /cOFIIF/ IF INT.AFINT.bFInT.CFInT.MMINT.TfINT 

c 

COMMON /INPUT/ I 0 * M A * * I C T A B . I PR I N T , I T W T AB • I * T A B . MZE T A . OX M A X , 

A EPSZ.FJ.G.GAMO.PO.PHlI .PIE.PRANOT.RBAK.S cALE.TO. 

K THETAI.TOLCFA.TOLZET.TOLZME.ZMUO.ZMVIS.ZNSTAN 

c 

COMMON /INTER/ CF aGT .CFaGP .CHPAR1 ,DX.D*RH0, HE .H* , 1BEG .HZETAM, 

A OOMZET.PhIP.PRE1Q3.RH0E.KH0UE.RMZETA,THETAP, 

K XI BASE .XIENO.ZETaTM.ZMZETA.ZMZETm.ZMZETP 

c 

COMMON /OUTPUT/ BDEL T A . CF , CH . DEL T A . DEL SO T , DELSTR , FL A T , F ORCE . H6 . 
A PE, PM I ,fi*,SUM<3DA,TE.THETA,T#.UE,X,XLARC,YR.Zl. 

K Z2,Z3.Z9,ZS.ZETA,ZME 

c 

COMMON /saved/ A.B.C,ZI1.ZI1P,ZI2,Z12P»Z13,ZI3P.ZI9.ZIS,ZI6,ZI7 
C 

ERASE1-PHI/THETA 
I F I N T» l 

DO 30 I » 1 , 50 

mmint»mzeta 

AFlNTM 

ZETA6-ZETA 

IF (ZETa . 6E • 1.0) GO TO 32 
BFINT ■ B 

cfint«c*zeta»zeta 

CALL INTZETIO., I..ZUPI 

BF I NT«B/ZETA 

CFInT»C 

CALL INTZETIO. , ZETA. ZI9) 

AF I N T» A *B 
BFlNT-O, 

call intzetizeTa»i..ziS> 

ZETA»(ErasEi/ZI IP* tZIR*ZISl )*«R mZETA 
60 TO 33 
32 BFINT»B/ZETA 
CFInT«C 

CALL INTZETIO., I. . Z I 1 ) 

BFlNT-8 

CF I nt«c*zeta»zEta 

ERASE2-1 ,/ZETA 

CALL INTZETIO. .ERASE2.Z12P) 

MMiNT»MZtTAM 
AF iNTmA.C 
CFINT-O. 

CALL INTZETIERASE2.1..ZI3P) 

ZETa«(ERaSE 1/|Z12P*ZI3P/ZETA)»ZI I )»«RMZETA 
33 DZETA ■ IZETA - ZETAGl/ZETAG 

IF lABS(OZETA) .LT. TOLZEt) 60 TO 3S 
IF I I . GE . 2) GO TO 76 
Z9« ZET A 
Z2-ZETAG 
GO TO 30 
76 Z3«Z‘< 


ZETA 1 

/COFI IF/ 

/INPUT/ 

/INPUT/ 

/INPUT/ 

/INTER/ 

/inter/ 

/INTER/ 

/OUTPUT/ 

/OUTPUT/ 

/OUTPUT/ 

/saved/ 

ZETA 16 
ZETA 17 
ZETA 18 
ZETA 19 
ZETA 20 
ZETA 21 


ZETA 29 
ZETA 2S 
ZETA 26 
ZETA 27 
ZETA 28 
ZETA 29 
ZETA 30 
ZETA 31 
ZETA 32 
ZETA 33 
ZETA 3H 
ZETA 3S 
ZETA 36 
ZETA 37 
ZETA 38 
ZETA 39 
ZETA 90 
ZETA 91 
ZETA 92 
ZETA 93 
ZETA 99 
ZETA 95 


ZETA 50 
ZETA 51 
ZETA 52 
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Zt»Z2 
Z8»ZET A 
Z2-ZETAG 

ZS*(Z8-Z3)/(22-Zl ) 
ZETA*(ZH-Z5»Z2>/( l.-ZS) 

30 CONTINUE 

»RITE(6,39) X. ZME, THETA, PHI 
38 FORMAT C S7H0»»ZETAIT FAILURE#.. 
1... / 22H0 AXIAL DISTANCE X >. 

2 E18.7. bX. 8HTHETAI •, EJH.7, bX, 


HR I TE ( 6 , bO ) Zl> Z2 » ZETa. Z3, Z8 
bO FORMAT ( 2QH ZETA (GUESSED) ■, 
1 ». 2EU.7 // > 

3b I F I NT « 2 

MM 1 NT »M2ET A 
AF I NT»A 
BFINT»B/ZETA 
CFINT-C 

ZETATM»ZETA»*ZMZETA 

IF (ZETA .GE. 1.0) GO TO 37 

CALL INTZET (0..ZETA.ZI6) 

AF I N T»A *b 
BF JnT.O. 

CALL INTZET(ZETa.I.,ZI 7) 
ERASE2»Zl8*ZIb 

DELS0T*(00MZET-Z16-ZI7)/ERASE2 
UELTA»THETA/ZMZET'A/ERASE2 
GO TO 38 

37 CALL INTZETIO.il. > Z I 2 ) 

MM I NT»MZETAM 
AFINT»A+C 
CF INT«0. 

CALL INTZETI 1. .ZETA.Z13) 
DELTA«THETA/ZMZETA/ZIl 
DELS0T*(ZETATM/ZMZETA-ZI3-ZI2)/ZI 1 
38 BDELTA » zETATM»OELTA 

delstr«theta»delsot 

return 

END 


ZETA S3 
ZETA 58 
ZETA SS 
ZETA 56 
ZETA 57 
ZETA B8 
ZETA S9 

SHAPE PARAMETER ITERATION FAlLUREZETA 60 
1PE18.7, bX, 11HMACH NO. - . ZETA 61 

6HPHII «, E1H.7 ) ZETA 62 

ZETA 63 

1 P3E 1 6 . 7 / 20h ZETA ( C AL C UL A TED ) ZE T A 68 

ZETA 6b 

ZETA 68 
ZETA 69 
ZETA 70 
ZETA 71 
ZETA 72 


ZETA 7b 
ZETA 76 
ZETA 77 
ZETA 78 
ZETA 79 
ZETA BO 
ZETA 81 
ZETA 82 
ZETA 83 
ZETA 88 
ZETA 8b 
ZETA 86 
ZETA 87 
ZETA 88 

ZETA 91 
ZETA 92 
ZETA 93 
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